/******************************************************************************************    AutoSphere:  A program to form a docking set of spheres from the 0 cluster **                 spheres genenerated by SPHGEN.****    Created: 11-8-99**    Author: Scott Pegg (Kuntz Group, Dept. of Pharm. Chem., UCSF)**            Jim Arnold (Kuntz Group, Dept. of Pharm. Chem., UCSF)****   Synopsis:****     This program reads in an input file given as a command line argument.**     This file specifies a SPHGEN output file, starting coordinates, an**     enclosed space, and output filenames (see the README for more info about**     the parameters) The spheres in "cluster 0" of the SPHGEN output file are**     selected with the following procedure:****     1. Remove spheres not within a box, large sphere, or close to a protein residue.**     2. Remove spheres distant from the protein surface.**     3. Remove spheres too close to other spheres, as such:**           x = sphere closest to the center of the enclosed space**           loop (until all spheres visited or removed):**              remove all spheres closer than given filt_cut to x**              x = next closest remaining sphere****     The spheres are treated as points--ie. the distance between two spheres**     is calculated as the distance between their centers. (I treat them this**     way because they're treated as points in DOCK.)****     Additions by Jim Arnold, Nov, 22, 1999.**        - Remove spheres distant from a protein residue.**        - Remove spheres distant from the protein surface.**        - Read molecule files using Dock4.0.1 functions (archived Dock).**        - Write out sphere sets in both .pdb and .sph formats.**        - Added debug and guards to the program.****************************************************************************************/#include  "AutoSphere.h"#include  <stdio.h>#include  <malloc.h>#include  <stdlib.h>#include  <string.h>#include  <math.h>#include  "logical.h"#include  "debug.h"#include  "/idk1/people/arnold/fulcrum/Dock4.0.1/source/dock/define.h"#include  "/idk1/people/arnold/fulcrum/Dock4.0.1/source/dock/global.h"#include  "/idk1/people/arnold/fulcrum/Dock4.0.1/source/dock/utility.h"#include  "/idk1/people/arnold/fulcrum/Dock4.0.1/source/dock/mol.h"#include  "/idk1/people/arnold/fulcrum/Dock4.0.1/source/dock/label.h"/**************** mark defs.**************/#define    DEBUG_FILE      "as.dbf"                            /** file with debug flags. **/#define    PARAM_DIR       "/idk1/people/arnold/fulcrum/Tsp"   /** directory for debug file. **//****************** mark globals****************/GLOBAL global = {0};            /* See global.h for global variables *//********************* mark prototypes*******************/void         read_args(FILE *file, TS_INFO *info);int          TestParams(TS_INFO *info);SPHERE*      read_spheres(TS_INFO *info, SPHERE *spheres);SPHERE*      AddLigSpheres(MOLECULE *mol, SPHERE *spheres);int          ReadMol(FILE_NAME molS, MOLECULE *mol_ref, MOLECULE *mol_init);void         SetLabel(LABEL *label);void         write_spheres(TS_INFO *info, SPHERE *spheres);SPHERE*      find_enclosed_spheres(TS_INFO *info, SPHERE *spheres, MOLECULE *mol);SPHERE*      FilterByResidue(TS_INFO *info, SPHERE *spheres, MOLECULE *mol);SPHERE*      SurfDistThin(TS_INFO *info, SPHERE *spheres, MOLECULE *mol);SPHERE*      greedy_thin(TS_INFO *info, SPHERE *spheres);SPHERE*      add_sphere(SPHERE *sphere, SPHERE *head);SPHERE*      remove_sphere(SPHERE *sphere, SPHERE *head);int          SphereGuard(SPHERE *spheres, STRING200 errorS);int          SphLen(SPHERE *sph);int          IsHeavy(ATOM *at);float        dist(SPHERE *a, SPHERE *b);void         show_info(TS_INFO *info);void         show_spheres(SPHERE *spheres);/** mark main **//****************************************************************************  Author:   Scott Pegg, Jim Arnold.**  Purpose:  Main routine for TS.**  Notes:    **************************************************************************/main(int argc, char *argv[]) {  const char   c[101];  FILE         *param_file;  TS_INFO      *info;  SPHERE       *spheres = NULL;  /* head of the list of spheres */  SPHERE       *temp, *temp2;  MOLECULE     mol_ref = {0}, mol_init = {0}, lig_ref = {0}, lig_init = {0};  if (argc == 1) {    printf("Usage: TS <parameter file>\n");    exit(0);  }  DebugInitI(PARAM_DIR, DEBUG_FILE);     /** read debug flags from file. **/  if (MONITOR) {printf("\nRunning AutoSphere.\n");}  param_file = fopen(argv[1], "r");  if (param_file == NULL) {    printf("Error: Could not open input file: %s\n", argv[1]);    exit(0);  }  info = malloc(sizeof(struct ts_info));  read_args(param_file, info);  if (TestParams(info) equal true)  {     spheres = read_spheres(info, spheres);     if (SphereGuard(spheres, "reading spheres, check the sph_file value in the .in file"))     {        if (MONITOR) {printf("  Read |%d| spheres from cluster 0 of file %s.\n", SphLen(spheres),info->sph_input);}        if (info->lig_sph_file[0] != NULL)	{           allocate_molecule(&lig_ref);     lig_ref.info.name = NULL;                allocate_molecule(&lig_init);    lig_init.info.name = NULL;              ReadMol(info->lig_sph_file, &lig_ref, &lig_init);                  if (MONITOR) {printf("  Read lig |%s| with %d atoms.\n",info->lig_sph_file, lig_ref.total.atoms);}           spheres = AddLigSpheres(&lig_ref, spheres);	}        allocate_molecule(&mol_ref);     mol_ref.info.name = NULL;             allocate_molecule(&mol_init);    mol_init.info.name = NULL;           ReadMol(info->receptor_filename, &mol_ref, &mol_init);               if (MONITOR) {printf("  Read receptor |%s| with %d atoms.\n",info->receptor_filename,mol_ref.total.atoms);}        spheres = find_enclosed_spheres(info, spheres, &mol_ref);        if (SphereGuard(spheres,             "finding spheres in volume, check residue or res_cut values in .in file - or remove TERs from .pdb"))        {           if (MONITOR) {printf("  Removed spheres outside defined volume:  %d remain.\n", SphLen(spheres));}           spheres = SurfDistThin(info, spheres, &mol_ref);           if (SphereGuard(spheres,                "surf dist filter, check surface_cutoff value in the .in file, or mol/sph coords"))  	   {              if (MONITOR){printf("  Removed spheres distant from receptor atoms:  %d remain.\n", SphLen(spheres));}               spheres = greedy_thin(info, spheres);              if (SphereGuard(spheres, "greedy filter, check the filter_cutoff value in the .in file"))   	      {                 if (MONITOR) {printf("  Filtered spheres using greedy algorithm:  %d remain.\n", SphLen(spheres));}                 if (MONITOR)                     {printf("  Writing %d spheres to .pdb and .sph: |%s|.\n",SphLen(spheres),info->sph_output);}                 write_spheres(info, spheres);  	      }	   }        }     }  }  free(info);  temp = spheres;  while (temp != NULL)   {    temp2 = temp->next;    free(temp);    temp = temp2;  }  if (MONITOR) {printf("  Exiting AutoSphere.\n\n");}}/** mark AddLigSpheres **//****************************************************************************  Author:   Jim Arnold**  Purpose:  Adds spheres from the ligand to the 0 cluster spheres.**  Notes:    **************************************************************************/SPHERE* AddLigSpheres(MOLECULE *mol, SPHERE *spheres){   int      i, ctI;   SPHERE   *new_sphere;   if (DB_LIG)       {printf("**AddLigSpheres, mol: %s, num sph: %d.\n", mol->info.name,SphLen(spheres));}   for (i = 0; i < mol->total.atoms; i++)   {      if (IsHeavy(&(mol->atom[i])))      {         new_sphere = malloc(sizeof(struct sphere));         new_sphere->x = mol->coord[i][0];          new_sphere->y = mol->coord[i][1];          new_sphere->z =mol->coord[i][2];         new_sphere->radius = 1.0;         new_sphere->visited = 0;         if (DB_LIG) {printf("  add lig sphere %d to the list.\n", ctI);}         spheres = add_sphere(new_sphere, spheres);         ctI++;     }   }   if (DB_LIG) {printf("  ret: %d sph from AddLigSpheres, added: %d.\n", SphLen(spheres), ctI);}   return (spheres);}/** mark SphereGuard **//****************************************************************************  Author:   Jim Arnold**  Purpose:  Generates an error message if the number of spheres is 0.**  Notes:    **************************************************************************/int SphereGuard(SPHERE *spheres, STRING200 errorS){   int retI;   if (DB_SPH_G) {printf("**SphereGuard, num sph: %d, error: |%s|.\n", SphLen(spheres), errorS);}   retI = SphLen(spheres);   if (retI equal 0)   {      printf("\n\nERROR:  %s.\n\n\n", errorS);   }   if (DB_SPH_G) {printf("  returning %d from SphereGuard.\n", retI);}   return (retI);}/** mark find_enclosed_spheres **//****************************************************************************  Author:   Scott Pegg**  Purpose:  A function to remove all spheres not within the specified **            enclosed area (a sphere or box).**  Notes:    **************************************************************************/SPHERE* find_enclosed_spheres(TS_INFO *info, SPHERE *spheres, MOLECULE *mol) {  SPHERE *temp, *temp2;  SPHERE center;  float x_lim1, x_lim2, y_lim1, y_lim2, z_lim1, z_lim2;  if (DB_ENCLOSED) {printf("**find_enclosed_spheres, num sph: %d.\n", SphLen(spheres));}  if (info->residue)  {     if (DB_ENCLOSED) {printf("  filter by residue num: %d.\n", info->res_num);}     spheres = FilterByResidue(info, spheres, mol);  }  else if (info->sphere)   {    if (DB_ENCLOSED) {printf("  filter by sphere.\n");}    /* make a dummy sphere for the center of the enclosure */    center.x = info->sph_x;    center.y = info->sph_y;    center.z = info->sph_z;    /* for each sphere, check if it's inside the enclosure radius, and        remove it if it's not */     temp = spheres;    while (temp != NULL)     {      temp2 = temp->next;      if (dist(&center, temp) > info->radius)       {	spheres = remove_sphere(temp, spheres);      }      temp = temp2;    }  }  else if (info->box)   {    if (DB_ENCLOSED) {printf("  filter with box.\n");}    /* set the x,y,z limits on the coordinates (Yes, I know I don't really       need to cast them explicity, but I like really clean compiles...) */    x_lim1 = info->sph_x + (float)(info->size_x / 2.0);    x_lim2 = info->sph_x - (float)(info->size_x / 2.0);    y_lim1 = info->sph_y + (float)(info->size_y / 2.0);    y_lim2 = info->sph_y - (float)(info->size_y / 2.0);    z_lim1 = info->sph_z + (float)(info->size_z / 2.0);    z_lim2 = info->sph_z - (float)(info->size_z / 2.0);    /* remove any spheres with a coordinate outside the limits */    temp = spheres;    while (temp != NULL)     {      temp2 = temp->next;       if (temp->x > x_lim1 || temp->x < x_lim2 || temp->y > y_lim1 ||	  temp->y < y_lim2 || temp->z > z_lim1 || temp->z < z_lim2)	spheres = remove_sphere(temp, spheres);      temp = temp2;    }  }  else  {     printf("ERROR:  No match for type of enclosed volume filtering.\n");  }  return spheres;}/** mark FilterByResidue **//****************************************************************************  Author:   Jim Arnold**  Purpose:  Removes spheres not within a given distance of a specified**            protein residue.**  Notes:    **************************************************************************/SPHERE* FilterByResidue(TS_INFO *info, SPHERE *spheres, MOLECULE *mol){   int             i, keep;   float           d;   SPHERE          *temp, *temp2;   static SPHERE   dsp;   if (DB_RES)    {      printf("**FilterByRes, res: %d, sph: %d nat: %d.\n",info->res_num,SphLen(spheres),mol->total.atoms);   }   if (SphLen(spheres) <= 0)   {      if (DB_RES) {printf("  GUARD:  %d spheres passed to FilterByResidue.\n", SphLen(spheres));}   }   else   {      temp = spheres;      while (temp != NULL)      {         keep = false;         for (i = 0; i < mol->total.atoms; i++)         {            if (DB_RES) {printf("  atom: %d, res num: %d.\n", i, mol->atom[i].res_num);}            if (mol->atom[i].res_num equal info->res_num)	    {               if (IsHeavy(&(mol->atom[i])))      	       {                  dsp.x = mol->coord[i][0]; dsp.y = mol->coord[i][1]; dsp.z = mol->coord[i][2];                  d = dist(temp, &dsp);                  if (DB_RES_L) {printf("    dist is %4.1f.\n", d);}                  if (d < info->res_cut)     	          {                     if (DB_RES_L) {printf("      dist < %4.1f, keep sph.\n", d);}                     keep = true;      	          }	       }            }	 }         if (keep equal false)         {            if (DB_RES_L) {printf("    keep is false, remove sphere.\n");}            temp2 = temp->next; /* temp get's freed if removed */            spheres = remove_sphere(temp, spheres);            temp = temp2;         }