voronoi.h File Reference

#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <ctype.h>
#include <sys/types.h>
#include <unistd.h>
#include "rpdb.h"
#include "writepdb.h"
#include "calc.h"
#include "utils.h"
#include "../src/qhull/qvoronoi.h"
#include "memhandler.h"

Go to the source code of this file.

Data Structures

struct  s_vvertice
struct  s_lst_vvertice

Defines

#define M_VORONOI_SUCCESS   0
#define M_APOLAR_AS   0
#define M_POLAR_AS   1
#define M_PREC_TOLERANCE   1e-5
#define M_BUFSIZE   1e7

Functions

s_lst_vverticeload_vvertices (s_pdb *pdb, int min_apol_neigh, float ashape_min_size, float ashape_max_size)
float testVvertice (float xyz[3], int curNbIdx[4], s_atm *atoms, float min_asph_size, float max_asph_size, s_lst_vvertice *lvvert)
void set_barycenter (s_vvertice *v)
int is_in_lst_vert (s_vvertice **lst_vert, int nb_vert, int v_id)
int is_in_lst_vert_p (s_vvertice **lst_vert, int nb_vert, s_vvertice *vert)
void write_pqr_vert (FILE *f, s_vvertice *v)
void write_pdb_vert (FILE *f, s_vvertice *v)
float get_verts_volume_ptr (s_vvertice **verts, int nvert, int niter, float correct)
void print_vvertices (FILE *f, s_lst_vvertice *lvvert)
void free_vert_lst (s_lst_vvertice *lvvert)


Define Documentation

#define M_APOLAR_AS   0

alpha sphere type - hydrophilic alpha sphere

Definition at line 57 of file voronoi.h.

Referenced by fill_vvertices(), get_vert_apolar_density(), set_descriptors(), updateIds(), write_pdb_vert(), and write_pqr_vert().

#define M_BUFSIZE   1e7

buffer size

Definition at line 61 of file voronoi.h.

#define M_POLAR_AS   1

tolerance for coordinate imprecion during alpha sphere search

Definition at line 58 of file voronoi.h.

Referenced by fill_vvertices().

#define M_PREC_TOLERANCE   1e-5

Definition at line 59 of file voronoi.h.

Referenced by testVvertice().

#define M_VORONOI_SUCCESS   0

alpha sphere type - hydrophobic alpha sphere

Definition at line 56 of file voronoi.h.

Referenced by load_vvertices().


Function Documentation

void free_vert_lst ( s_lst_vvertice lvvert  ) 

## FUNCTION: free_vert_lst

## SPECIFICATION: Free memory

## PARAMETERS: @ s_lst_vvertice *lvvert : Data to free

## RETURN: void

Definition at line 530 of file voronoi.c.

References s_lst_vvertice::h_tr, my_free(), s_lst_vvertice::pvertices, s_lst_vvertice::tr, and s_lst_vvertice::vertices.

Referenced by c_lst_pocket_free().

00531 {
00532     if(lvvert) {
00533 
00534         if(lvvert->vertices) {
00535             my_free(lvvert->vertices) ;
00536             lvvert->vertices = NULL ;
00537         }
00538         if(lvvert->pvertices) {
00539             my_free(lvvert->pvertices) ;
00540             lvvert->pvertices = NULL ;
00541         }
00542         if(lvvert->tr) {
00543                 my_free(lvvert->tr) ;
00544                 lvvert->tr = NULL ;
00545         }
00546         if(lvvert->h_tr) {
00547                 my_free(lvvert->h_tr) ;
00548                 lvvert->h_tr = NULL ;
00549         }
00550         my_free(lvvert) ;
00551     }
00552 }

float get_verts_volume_ptr ( s_vvertice **  verts,
int  nvert,
int  niter,
float  correct 
)

## FUNCTION: get_verts_volume_ptr

## SPECIFICATION: Get an monte carlo approximation of the volume occupied by the alpha spheres given in argument (list of pointers)

## PARAMETRES: @ s_vvertice **verts: List of pointer to alpha spheres @ int nvert: Number of spheres @ int niter: Number of monte carlo iteration to perform @ float correct: radius for which the size of an alpha sphere should be corrected in order to calculate the volume

## RETURN: float: volume.

Definition at line 454 of file voronoi.c.

References rand_uniform(), s_vvertice::ray, s_vvertice::x, s_vvertice::y, and s_vvertice::z.

Referenced by set_descriptors().

00455 {
00456         int i = 0, j = 0,
00457                 nb_in = 0;
00458 
00459         float xmin = 0.0, xmax = 0.0,
00460                   ymin = 0.0, ymax = 0.0,
00461                   zmin = 0.0, zmax = 0.0,
00462                   xtmp = 0.0, ytmp = 0.0, ztmp = 0.0,
00463                   xr   = 0.0, yr   = 0.0, zr   = 0.0,
00464                   vbox = 0.0 ;
00465 
00466         s_vvertice *vcur = NULL ;
00467 
00468         /* First, search extrems coordinates to get a contour box of the molecule */
00469         for(i = 0 ; i < nvert ; i++) {
00470                 vcur = verts[i] ;
00471 
00472                 if(i == 0) {
00473                         xmin = vcur->x - vcur->ray + correct ; xmax = vcur->x + vcur->ray + correct ;
00474                         ymin = vcur->y - vcur->ray + correct ; ymax = vcur->y + vcur->ray + correct ;
00475                         zmin = vcur->z - vcur->ray + correct ; zmax = vcur->z + vcur->ray + correct ;
00476                 }
00477                 else {
00478                 /* Update the minimum and maximum extreme point */
00479                         if(xmin > (xtmp = vcur->x - vcur->ray + correct)) xmin = xtmp ;
00480                         else if(xmax < (xtmp = vcur->x + vcur->ray + correct)) xmax = xtmp ;
00481 
00482                         if(ymin > (ytmp = vcur->y - vcur->ray + correct)) ymin = ytmp ;
00483                         else if(ymax < (ytmp = vcur->y + vcur->ray + correct)) ymax = ytmp ;
00484 
00485                         if(zmin > (ztmp = vcur->z - vcur->ray + correct)) zmin = ztmp ;
00486                         else if(zmax < (ztmp = vcur->z + vcur->ray + correct)) zmax = ztmp ;
00487                 }
00488         }
00489 
00490         /* Next calculate the contour box volume */
00491         vbox = (xmax - xmin)*(ymax - ymin)*(zmax - zmin) ;
00492 
00493         /* Then apply monte carlo approximation of the volume.   */
00494         for(i = 0 ; i < niter ; i++) {
00495                 xr = rand_uniform(xmin, xmax) ;
00496                 yr = rand_uniform(ymin, ymax) ;
00497                 zr = rand_uniform(zmin, zmax) ;
00498 
00499                 for(j = 0 ; j < nvert ; j++) {
00500                         vcur = verts[j] ;
00501                         xtmp = vcur->x - xr ;
00502                         ytmp = vcur->y - yr ;
00503                         ztmp = vcur->z - zr ;
00504 
00505                 /* Compare r^2 and dist(center, random_point)^2 */
00506                         if(((correct+vcur->ray)*(vcur->ray + correct)) > (xtmp*xtmp + ytmp*ytmp + ztmp*ztmp)) {
00507                         /* The point is inside one of the vertice!! */
00508                                 nb_in ++ ; break ;
00509                         }
00510                 }
00511         }
00512 
00513         /* Ok lets just return the volume Vpok = Nb_in/Niter*Vbox */
00514         return ((float)nb_in)/((float)niter)*vbox;
00515 }

int is_in_lst_vert ( s_vvertice **  lst_vert,
int  nb_vert,
int  v_id 
)

## FUNCTION: is_in_lst_vert

## SPECIFICATION: Says wether a vertice of id v_id is in a list of vertices or not

## PARAMETRES:

## RETURN: 1 if the vertice is in the tab, 0 if not

Definition at line 567 of file voronoi.c.

Referenced by get_mol_vert_neigh().

00568 {
00569         int i ;
00570         for(i = 0 ; i < nb_vert ; i++) {
00571                 if(v_id == lst_vert[i]->id) return 1 ;
00572         }
00573 
00574         return 0 ;
00575 }

int is_in_lst_vert_p ( s_vvertice **  lst_vert,
int  nb_vert,
s_vvertice vert 
)

## FUNCTION: is_in_lst_vert

## SPECIFICATION: Says wether a vertice of id v_id is in a list of vertices or not

## PARAMETRES:

## RETURN: 1 if the vertice is in the tab, 0 if not

Definition at line 590 of file voronoi.c.

00591 {
00592         int i ;
00593         for(i = 0 ; i < nb_vert ; i++) {
00594                 if(vert == lst_vert[i]) return 1 ;
00595         }
00596 
00597         return 0 ;
00598 }

s_lst_vvertice* load_vvertices ( s_pdb pdb,
int  min_apol_neigh,
float  asph_min_size,
float  asph_max_size 
)

## FUNCTION: s_lst_vvertice

## SPECIFICATION: Calculate voronoi vertices using an ensemble of atoms, and then load resulting vertices into a s_lst_vvertice structure. The function call an external programm qvoronoi, part of qhull programme which can be download at: http://www.qhull.org/download/ or installed with apt-get install qhull-bin

## PARAMETRES: @ s_pdb *pdb : PDB informations @ int min_apol_neigh : Number of apolar neighbor of a vertice to be considered as apolar @ float asph_min_size : Minimum size of voronoi vertices to retain @ float asph_max_size : Maximum size of voronoi vertices to retain

## RETURN: s_lst_vvertice * :The structure containing the list of vertices.

Definition at line 90 of file voronoi.c.

References fill_vvertices(), s_lst_vvertice::h_tr, s_pdb::latoms, M_VORONOI_SUCCESS, my_free(), my_malloc(), my_realloc(), s_lst_vvertice::n_h_tr, s_pdb::natoms, s_atm::symbol, s_atm::x, s_atm::y, and s_atm::z.

Referenced by search_pocket().

00091 {
00092         int i, nb_h=0;
00093         s_atm *ca = NULL ;
00094         s_lst_vvertice *lvvert = NULL ;
00095         char tmpn1[250]="";
00096         char tmpn2[250]="";
00097         pid_t pid=getpid();
00098         sprintf(tmpn1,"/tmp/qvoro_tmp_%d.dat",pid);
00099         sprintf(tmpn2,"/tmp/fpocket_qvor_%d.dat",pid);
00100         FILE *fvoro = fopen(tmpn1, "w+");
00101         FILE *ftmp=fopen(tmpn2,"w");
00102         /*FILE *ftmp=fopen("/tmp/fpocket_qvor.dat","w");
00103         FILE *fvoro = fopen("/tmp/qvoro_tmp.dat", "w+");*/
00104 /*      lvvert->h_tr=(int *)my_malloc(sizeof(int));*/
00105 
00106 
00107 
00108         if(fvoro != NULL) {
00109                 lvvert = (s_lst_vvertice *)my_malloc(sizeof(s_lst_vvertice)) ;
00110                 lvvert->h_tr=NULL;
00111                 /* Loop a first time to get out how many heavy atoms are in the file */
00112                 for(i = 0; i <  pdb->natoms ; i++){
00113                         ca = (pdb->latoms)+i ;
00114                         if(strcmp(ca->symbol,"H")) {
00115                                 lvvert->h_tr=(int *)my_realloc(lvvert->h_tr,sizeof(int)*(i-nb_h+1)) ;
00116                                 lvvert->h_tr[i-nb_h]=i ;
00117                         }
00118                         else nb_h++;
00119                 }
00120                 lvvert->n_h_tr=i-nb_h;
00121 
00122                 /* Write the header for qvoronoi */
00123                 fprintf(fvoro,"3 rbox D3\n%d\n", lvvert->n_h_tr);
00124                 /* Loop a second time for the qvoronoi input coordinates */
00125                 for(i = 0; i <  pdb->natoms ; i++){
00126                         ca = (pdb->latoms)+i ;
00127                         if(strcmp(ca->symbol,"H")) {
00128                         /* Only if this is a heavy atom export it for voronoi tesselation,
00129                          * else discard it */
00130                                 fprintf(fvoro,"%.3f %.3f %.3f \n", ca->x, ca->y, ca->z);
00131                         }
00132                 }
00133 
00134                 fflush(fvoro) ;
00135                 rewind(fvoro);
00136                 
00137                 //int status = system("qvoronoi p i Pp Fn < voro_tmp.dat > voro.tmp") ;
00138                 run_qvoronoi(fvoro,ftmp);
00139                 int status=M_VORONOI_SUCCESS;
00140 
00141 
00142               
00143                 if(status == M_VORONOI_SUCCESS) {
00144                         fill_vvertices(lvvert, tmpn2, pdb->latoms, pdb->natoms,
00145                                                         min_apol_neigh, asph_min_size, asph_max_size);
00146                 }
00147                 else {
00148                         my_free(lvvert);
00149                         lvvert = NULL ;
00150                         fprintf(stderr, "! Voronoi command failed with status %d...\n", status) ;
00151                 }
00152                 
00153           
00154         }
00155         else {
00156                 fprintf(stderr, "! File for Voronoi vertices calculation couldn't be opened...\n") ;
00157         }
00158         fclose(fvoro) ;
00159         fclose(ftmp);
00160         remove(tmpn1);
00161         remove(tmpn2);
00162         return lvvert ;
00163 }

void print_vvertices ( FILE *  f,
s_lst_vvertice lvvert 
)

## FUNCTION: print_vvertices

## SPECIFICATION: Print function.

## PARAMETERS: @ FILE *f : Buffer to print in @ s_lst_vvertice *lvvert : Vertices to print

## RETURN: void

Definition at line 411 of file voronoi.c.

References dist(), s_atm::id, s_vvertice::neigh, s_lst_vvertice::nvert, s_vvertice::sort_x, s_lst_vvertice::vertices, s_atm::x, s_vvertice::x, s_atm::y, s_vvertice::y, s_atm::z, and s_vvertice::z.

00412 {
00413         if(lvvert) {
00414                 if(lvvert->vertices) {
00415                         int i ;
00416                         for(i = 0 ; i < lvvert->nvert ; i++) {
00417                                 s_vvertice *v = &(lvvert->vertices[i]) ;
00418                                 if( v->neigh[0] &&  v->neigh[1] && v->neigh[2] &&  v->neigh[3]) {
00419                                         fprintf(f, "====== Vertice %d: =====\n", i);
00420                                         fprintf(f, "- x = %f, y = %f, z = %f\n", v->x, v->y, v->z);
00421                                         fprintf(f, "- ix = %d\n",v->sort_x);
00422 
00423                                         float d1 = dist(v->x, v->y, v->z, v->neigh[0]->x, v->neigh[0]->y, v->neigh[0]->z) ;
00424                                         float d2 = dist(v->x, v->y, v->z, v->neigh[1]->x, v->neigh[1]->y, v->neigh[1]->z) ;
00425                                         float d3 = dist(v->x, v->y, v->z, v->neigh[2]->x, v->neigh[2]->y, v->neigh[2]->z) ;
00426                                         float d4 = dist(v->x, v->y, v->z, v->neigh[3]->x, v->neigh[3]->y, v->neigh[3]->z) ;
00427 
00428                                         fprintf(f, "- Neighbour: \n1 - %f (%f %f %f: %d) \n2 - %f (%f %f %f: %d)\n3 - %f (%f %f %f: %d)\n4 - %f (%f %f %f: %d)\n", d1, v->neigh[0]->x, v->neigh[0]->y, v->neigh[0]->z, v->neigh[0]->id, d2, v->neigh[1]->x, v->neigh[1]->y, v->neigh[1]->z, v->neigh[1]->id, d3, v->neigh[2]->x,  v->neigh[2]->y ,  v->neigh[2]->z, v->neigh[2]->id, d4, v->neigh[3]->x, v->neigh[3]->y, v->neigh[3]->z, v->neigh[3]->id);
00429                                 }
00430                         }
00431                 }
00432         }
00433 }

void set_barycenter ( s_vvertice v  ) 

## FUNCTION: set_barycenter

## SPECIFICATION: Set barycenter of a vertice using it's 4 contacting atoms.

## PARAMETERS: @ s_vvertice *v: The vertice

## RETURN: void

Definition at line 321 of file voronoi.c.

References s_vvertice::bary, s_vvertice::neigh, s_atm::x, s_atm::y, and s_atm::z.

Referenced by fill_vvertices().

00322 {
00323         int i ;
00324         float xsum = 0.0,
00325                   ysum = 0.0,
00326                   zsum = 0.0 ;
00327 
00328         for(i = 0 ; i < 4 ; i++) {
00329                 xsum += v->neigh[i]->x ;
00330                 ysum += v->neigh[i]->y ;
00331                 zsum += v->neigh[i]->z ;
00332         }
00333 
00334         v->bary[0] = xsum*0.25 ;
00335         v->bary[1] = ysum*0.25 ;
00336         v->bary[2] = zsum*0.25 ;
00337 
00338 }

float testVvertice ( float  xyz[3],
int  curNbIdx[4],
s_atm atoms,
float  min_asph_size,
float  max_asph_size,
s_lst_vvertice lvvert 
)

## FUNCTION: testVvertice

## SPECIFICATION: Test if alpha sphere conditions are fulfilled for current vertice

## PARAMETERS: @ float xyz[3] : Coordinates of current vertice @ int curNbIdx[4] : Indexes of atomic neighbours of the current vertice @ s_atm *atoms : List of all atoms @ float min_asph_size : Minimum size of alpha spheres. @ float max_asph_size : Maximum size of alpha spheres.

## RETURN: float : -1 if conditions are not fulfilled, else the alpha sphere radius is returned.

Definition at line 359 of file voronoi.c.

References dist(), s_lst_vvertice::h_tr, M_PREC_TOLERANCE, s_atm::x, s_atm::y, and s_atm::z.

Referenced by fill_vvertices().

00362 {
00363         float x = xyz[0],
00364                   y = xyz[1],
00365                   z = xyz[2] ;
00366 
00367         s_atm *cura = &(atoms[lvvert->h_tr[curNbIdx[0]]]) ;
00368 
00369         float distVatom1 = dist(x, y, z, cura->x, cura->y, cura->z) ;
00370         float distVatom2,
00371                   distVatom3,
00372                   distVatom4;
00373 
00374         if(min_asph_size <= distVatom1 + M_PREC_TOLERANCE  && distVatom1 - M_PREC_TOLERANCE <= max_asph_size){
00375                 cura = &(atoms[lvvert->h_tr[curNbIdx[1]]]) ;
00376                 distVatom2 = dist(x, y, z, cura->x, cura->y, cura->z);
00377 
00378                 cura = &(atoms[lvvert->h_tr[curNbIdx[2]]]) ;
00379                 distVatom3 = dist(x, y, z, cura->x, cura->y, cura->z);
00380 
00381                 cura = &(atoms[lvvert->h_tr[curNbIdx[3]]]) ;
00382                 distVatom4=dist(x, y, z, cura->x, cura->y, cura->z);
00383 
00384                 /* Test if all 4 neighbours are on the alpha sphere surface
00385                  * (approximate test) */
00386                 if(fabs(distVatom1-distVatom2) < M_PREC_TOLERANCE &&
00387                    fabs(distVatom1-distVatom3)  < M_PREC_TOLERANCE &&
00388                    fabs(distVatom1-distVatom4) < M_PREC_TOLERANCE){
00389 
00390                         return distVatom1;
00391                 }
00392 
00393         }
00394         return -1.0;
00395 }

void write_pdb_vert ( FILE *  f,
s_vvertice v 
)

-- --

OUTPUT FUNCTIONS

-- -- -- ## FUNCTION: void write_pdb_vertice(FILE *f, s_vvertice *v)

## SPECIFICATION: Write a voronoi vertice in pdb format.

## PARAMETRES: @ FILE *f: file to write in @ s_vvertice *v: The vertice

## RETURN:

Definition at line 624 of file voronoi.c.

References s_vvertice::id, M_APOLAR_AS, s_vvertice::resid, s_vvertice::type, write_pdb_atom_line(), s_vvertice::x, s_vvertice::y, and s_vvertice::z.

Referenced by write_pockets_single_pdb().

00625 {
00626         if(v->type==M_APOLAR_AS) write_pdb_atom_line(f, "HETATM", v->id, "APOL",
00627                                                                                                  ' ', "STP", "C", v->resid, ' ',
00628                                                                                                  v->x, v->y, v->z, 0.0, 0.0,
00629                                                                                                  "Ve", -1);
00630 
00631         else write_pdb_atom_line(f, "HETATM", v->id, " POL", ' ', "STP", "C",
00632                                                                 v->resid, ' ',v->x, v->y, v->z,0.0, 0.0,
00633                                                                 "Ve", -1);
00634 }

void write_pqr_vert ( FILE *  f,
s_vvertice v 
)

## FUNCTION: write_pqr_vertice

## SPECIFICATION: Write a voronoi vertice in pqr format.

## PARAMETRES: @ FILE *f : file to write in @ s_vvertice *v : The vertice

## RETURN: void

Definition at line 651 of file voronoi.c.

References s_vvertice::id, M_APOLAR_AS, s_vvertice::ray, s_vvertice::resid, s_vvertice::type, write_pqr_atom_line(), s_vvertice::x, s_vvertice::y, and s_vvertice::z.

Referenced by mdpocket_characterize(), write_mdpockets_concat_pqr(), write_pocket_pqr(), write_pocket_pqr_DB(), and write_pockets_single_pqr().

00652 {
00653         if(v->type==M_APOLAR_AS) write_pqr_atom_line(f, "ATOM", v->id, "APOL", ' ',
00654                                                                                                  "STP", " ", v->resid, ' ',
00655                                                                                                   v->x, v->y, v->z, 0.0, v->ray);
00656 
00657         else write_pqr_atom_line(f, "ATOM", v->id, " POL", ' ', "STP", " ",
00658                                                          v->resid, ' ',v->x, v->y, v->z,0.0, v->ray);
00659 }


Generated on Mon Jun 7 16:44:23 2010 for fpocket by  doxygen 1.5.6