automolprep.perl

最经典的分子对结软件

   if ($DB_AUTO_SPH_WR) {printf("**WriteAutoSphereInput, file: $asph_in\n");}

   open(asphF, ">$a_sph_in") || die "cannot open $a_sph_in for writing.";
   if ($DB_AUTO_SPH_WR) {printf("  opened sph_in file: $a_sph_in\n");}

   $res_cut = "10.0";
   $fil_cut = "1.4";
   $surf_cut = "6.0";
   $sphere = "0";
   $rad = "1.0";
   $box = "0";
   $box_center = "0";
   $rec_file = sprintf("%s%s%s", $root_file_name, $extn[0], ".pdb");
   $sph_file = sprintf("%s%s%s", $root_file_name, $extn[0], $sph_ex);
   $out_file = sprintf("%s%s%s", $root_file_name, $extn[0], "_as");

   printf asphF ("residue         1 \n");
   printf asphF ("res_num         $cent_res \n");
   printf asphF ("residue_cutoff  $res_cut \n");
   printf asphF ("sphere          $sphere \n");
   printf asphF ("radius          $rad \n");
   printf asphF ("box             $box \n");
   printf asphF ("box_center      $box_center \n");
   printf asphF ("filter_cutoff   $fil_cut \n");
   printf asphF ("surface_cutoff  $surf_cut \n");
   printf asphF ("receptor_file   $rec_file \n");
   printf asphF ("sph_file        $sph_file \n");
   printf asphF ("output_file     $out_file \n");

   close(asphF);

   if ($DB_AUTO_SPH_WR) {printf("  returning from WriteAutoSphereInput, wrote: $a_sph_in.\n");}
}

# mark CallSPHGEN
##############################################################################
# Author:   Jim Arnold
# Purpose:  Calls SPHGEN to generate spheres for a given receptor molecule.
# Notes:    
##############################################################################
sub CallSPHGEN($root_file_name, $surf_ex, $sph_ex)
{
   if ($DB_SPHGEN) {printf("**CallSPHGEN for file: $root_file_name\n");}

   printf("Deleting any left over sphgen temp files.\n");
   system("rm OUTSPH");
   system("rm temp1.*");   
   system("rm temp2.*");   
   system("rm temp3.*");
   
   $sph_in = "INSPH";

   WriteSPHGEN_Input($root_file_name, $sph_in, $surf_ex, $sph_ex);

   if (system("sphgen $sph_in") != 0)
   {
      printf("\n--> ERROR:  system call to SPHGEN with input $sph_in was NOT successful.\n\n");
   }

   unlink("$sph_in");
   unlink("OUTSPH");
   if ($DB_SPHGEN) {printf("  returning from having called SPHGEN.\n");}
}

# mark WriteSPHGEN_Input
##############################################################################
# Author:   Jim Arnold
# Purpose:  Writes the SPHGEN input file.
# Notes:    Params can later be set from a preferences file.
##############################################################################
sub WriteSPHGEN_Input($root_file_name, $sph_in, $surf_ex, $sph_ex)
{
   if ($DB_SPHGEN_WR) {printf("**WriteSPHGEN_Input for file: $root_file_name\n");}

   open(sphF, ">$sph_in") || die "cannot open $sph_in for writing.";
   if ($DB_SPHGEN_WR) {printf("  opened sph_in file: $sph_in\n");}

   $surf_file = sprintf("%s%s%s", $root_file_name, $extn[0], $surf_ex);
   $srftp = "R";
   $dentag = "X";
   $dotlim = "0.0";
   $radmax = "4.0";
   $radmin = "1.4";
   $sph_file = sprintf("%s%s%s", $root_file_name, $extn[0], $sph_ex);
   
   unlink($sph_file);

   printf sphF ("$surf_file \n");
   printf sphF ("$srftp \n");
   printf sphF ("$dentag \n");
   printf sphF ("$dotlim \n");
   printf sphF ("$radmax \n");
   printf sphF ("$radmin \n");
   printf sphF ("$sph_file \n");

   close(sphF);

   if ($DB_SPHGEN_WR) {printf("  returning from having called WriteSPHGEN_Input.\n");}
}

# mark MakeSurface
##############################################################################
# Author:   Jim Arnold
# Purpose:  Creates a molecular surface upon the receptor site molecule by
#           calling the designated program.
# Notes:    
##############################################################################
sub MakeSurface($root_file_name, $surf_ex)
{
   if ($DB_SURF) {printf("**MakeSurface for mol: $root_file_name\n");}

   $dir = "/idk2/private/midas/bin";
   $prog = "dms";
   $opt = " -n -a ";

   $surf_file = sprintf("%s%s%s", $root_file_name, $extn[0], $surf_ex);
   $mol_name = sprintf("%s%s%s", $root_file_name, $extn[0], ".pdb");

   $command = sprintf("%s%s%s%s%s%s%s%s", $dir, "/", $prog, " ", $mol_name, $opt, "-o ", $surf_file);

   if ($DB_SURF) {printf("  giving dms command to system: %s.\n", $command);}
   system("$command");

   if ($DB_SURF) {printf("  returning from making surface.\n");}
}

# mark PrepMols
##################################################################################
# Author:   Jim Arnold
# Purpose:  Writes a Sybyl .spl script and uses it to place hydrogens and
#           charges on the given molecule.
# Notes:    The read_pdb function in Dock4.0.1 assumes that the protein .pdb
#           file is not interrupted by "TER" fields.  If the atoms in the
#           protein are interrupted, read_pdb will only read to the first
#           "TER".  So, the .pdb written out by PrepMols has the TER's
#           removed so the full protein can be processed by Dock and AutoSphere.
##################################################################################
sub PrepMols($root_file_name)
{
   if ($DB_SYB_PREP) {printf("**SybylPrep, file: %s.\n", $root_file_name);}

   $syb_nameS = "temp.spl";

   $e = 0;
   while ($e < 3)
   {
      $ext = $extn[$e];
      if ($DB_SYB_PREP) {printf("  Sybyl prep component: $ext.\n");}

      WriteSPL($syb_name, $root_file_name, $ext);
      if (system("/idk2/private/sybyl/trigo/trigo sybyl6.5 < temp.spl > syb_out") != 0)
      {
         printf("\n--> ERROR:  system call to Sybyl from SybylPrep was NOT successful.\n\n");
      }
      unlink("syb_out");

      if ($ext eq "_rec")
      {
         $pdb_name = sprintf("%s%s%s", $root_file_name, $ext, ".pdb");
         $key = "TER";
         $tfile = "temp_pdb_file.pdb";
         RemoveKeyFromPDB($pdb_name, $key, $tfile);

         if ($DB_SYB_PREP) {printf("  moving $tfile to $pdb_name.\n");}
         system("mv $tfile $pdb_name");
      }

      $e++;
   }

   unlink($syb_nameS);
   if ($DB_SYB_PREP) {printf("  returning from SybylPrep.\n");}
}

# mark RemoveKeyFromPDB
##############################################################################
# Author:   Jim Arnold
# Purpose:  Removes lines that begin with the passed Key.
# Notes:    
##############################################################################
sub RemoveKeyFromPDB($pdb_name, $key, $tfile)
{
   if ($DB_RM_KEY) {printf("**RemoveKeyFromPDB, file: |$pdb_name|, key: |$key|, tfile: |$tfile|.\n");}

   open(t_inF, "$pdb_name") || die "cannot open $pdb_name for read.";
   if ($DB_RM_KEY) {printf("  opened pdb file: |%s| for read.\n", $pdb_name);}
   open(t_outF, ">$tfile") || die "cannot open $tfile for write.";
   if ($DB_RM_KEY) {printf("  opened pdb file: |%s| for write.\n", $tfile);}

   while(<t_inF>)
   {
      ($s1, $s2) = split(/\s+/,$_);

      if ($s1 eq $key)
      {
         if ($DB_RM_KEY) {printf("  Not writing line: |$_|.\n");}
      }
      else
      {
         printf t_outF ("%s", $_);
      }
   }
   close(t_inF);
   close(t_outF);

   if ($DB_RM_KEY) {printf("  returning from RemoveKeyFromPDB.\n");}
}

# mark WriteSPL
##############################################################################
# Author:   Jim Arnold
# Purpose:  Writes a Sybyl .spl file that is used to read in .pdb files, put
#           hydrogens and charges on them, and write out .mol2 files.
# Notes:    $pct = "%";
##############################################################################
sub WriteSPL($syb_name, $root_file_name, $ext)
{
   if ($DB_WRITE_SPL) {printf("**WriteSPL, root file: %s, ext: %s.\n", $root_file_name, $ext);}

   unlink($syb_nameS);
   open(outF, ">$syb_nameS") || die "cannot open syb_nameS for writing.";
   if ($DB_WRITE_SPL) {printf("  opened sybyl .spl file: |%s|.\n", $syb_nameS);}

   $mol_pdb = sprintf("%s%s%s", $root_file_name, $ext, ".pdb");
   $mol_mol2 = sprintf("%s%s%s", $root_file_name, $ext, ".mol2");

   printf outF ("echo temp.spl being executed for $mol_pdb and $mol_mol2\n");
   printf outF ("   brookhaven in m1 $mol_pdb NO \n");
   printf outF ("   brookhaven out m1 $mol_pdb \n");
   printf outF ("   command BIO ADDH m1 ALL \n");
   printf outF ("   biopolymer load charges m1 KOLL_ALL \n");
   printf outF ("   mol2 out m1 $mol_mol2 \n");

   close(outF);
   if ($DB_WRITE_SPL) {printf("  returning from WriteSPL.\n");}
}

# mark ParsePDB
##############################################################################
# Author:   Jim Arnold
# Purpose:  Parses the .pdb file into receptor, ligand and waters.
# Notes:    
##############################################################################
sub ParsePDB($pdb_file)
{
   if ($DB_PARSE) {printf("**ParsePDB, pdb: |$pdb_file|.\n");}

   if ($DB_PARSE) {printf("  Opening pdb file: |%s|.\n\n",$pdb_file);}
   open(pdbF, "$pdb_file") || die "cannot open $pdb_file for reading.";

   $ra = 0; $la = 0; $wa = 0; $ma = 0;
   $rb = 0; $lb = 0; $wb = 0; $mb = 0;

   while(<pdbF>)
   {
      chop($_);
      ($s1, $s2, $s3, $s4, $s5) = split(/\s+/,$_);

      if ($s1 eq "HEADER")
      {
         if ($DB_PARSE) {printf("    header is: |%s|.\n", $_);}
         $header = $_;
      }
      elsif ($s1 eq "ATOM")
      {
         if ($s5 eq "I")
         {
            if ($DB_PARSE_LIG) {printf("    have ligand atom: %s %s %s %s.\n", $s2, $s3, $s4, $s5);}
            $l[$la] = sprintf("%s", $_);
            if ($DB_PARSE_LIG) {printf("    count: %d, stored line: |%s|.\n", $la, $l[$la]);}
            $la++;
         }
         else
         {
            if ($DB_PARSE_PROT) {printf("    have receptor atom: %s %s %s %s.\n", $s2, $s3, $s4, $s5);}
            $r[$ra] = $_;
            $ra++;
         }
      }
      elsif ($s1 eq "HETATM")
      {
         if (($s4 eq "HOH") || ($s4 eq "WAT"))
         {
            if ($DB_PARSE_WAT) {printf("    have water: %s %s %s %s.\n", $s2, $s3, $s4, $s5);}
            $w[$wa] = $_;
            $wa++;
         }
         elsif (($s3 eq "Ca") || ($s3 eq "MG") || ($s3 eq "Mg") ||
                ($s3 eq "NA") || ($s3 eq "Na") || ($s3 eq "ZN") || ($s3 eq "Zn") ||
                ($s3 eq "FE") || ($s3 eq "Fe") ||  ($s3 eq "K"))
         {
            if ($DB_PARSE_PROT) {printf("    have receptor metal atom: %s %s %s %s.\n", $s2, $s3, $s4, $s5);}
            $r[$ra] = $_;
            $ra++;
         }
         else
         {
            if ($DB_PARSE_LIG) {printf("    have ligand atom: %s %s %s %s.\n", $s2, $s3, $s4, $s5);}
            $l[$la] = $_;
            $la++;
         }
      }
      else
      {
         ## skip this line.
      }
   }

   if ($DB_PARSE) {printf("  Num Atoms, rec: %d, lig: %d, wat: %d, met: %d.\n", $ra, $la, $wa, $ma);}

   if ($ra > 0)
   {
      $ext = $extn[0]; 
      $max = $ra; 
      $i = 0;
      while ($i < $max)
         {$lines[$i] = $r[$i]; $i++;}
      WritePDB($pdb_file, $ext, $lines, $max, $header);
   }

   if ($la > 0)
   {
      $ext = $extn[1]; 
      $max = $la; 
      $i = 0;
      while ($i < $max)
         {$lines[$i] = $l[$i]; $i++;}
      WritePDB($pdb_file, $ext, $lines, $max, $header);
   }

   if ($wa > 0)
   {
      $ext = $extn[2]; 
      $max = $wa; 
      $i = 0;
      while ($i < $max)
         {$lines[$i] = $w[$i]; $i++;}
      WritePDB($pdb_file, $ext, $lines, $max, $header);
   }

   close (pdbF);
   if ($DB_PARSE) {printf("  Returning from ParsePDB file.\n");}
}

# mark WritePDB
##############################################################################
# Author:   Jim Arnold
# Purpose:  Writes a .pdb file.
# Notes:    
##############################################################################
sub WritePDB($pdb_file, $ext, $lines, $max, $header)
{
   if ($DB_WRITE) {printf("**WritePDB, file: |%s|, exten: |%s|, max: %d.\n", $pdb_file, $ext, $max);}
   $BIG_NUM = 100000;
   
   ($root) = split(/.pdb/,$pdb_file);
   $outS = sprintf("%s%s%s", $root, $ext, ".pdb");
   open(outF, ">$outS") || die "cannot open outS for writing.";
   if ($DB_WRITE) {printf("  opened pdb output file: |%s|.\n", $outS);}

   printf outF ("REMARK   %s created by AutoMolPrep.\n", $outS);

   $i = 0;
   while ($i < $max)
   {
      if ($DB_WRITE) {printf("  parse line %d is: |%s|.\n", $i, $lines[$i]);}

      ($s1,$s2,$s3) = split(/\s+/,$lines[$i]);

      if (length($s3) > 4)
      {
         ($s1,$s2,$s_t,$s5,$s6,$s7,$s8,$s9,$s10,$s11) = split(/\s+/,$lines[$i]);

         $s3 = substr($s_t, 0, 4);
         $s4 = substr($s_t, 4);
         if ($DB_WRITE_B) {printf("  s3: |$s3|, s4: |$s4|.\n");}
      }
      else
      {
         ($s1,$s2,$s3,$s4,$s5,$s6,$s7,$s8,$s9,$s10,$s11) = split(/\s+/,$lines[$i]);
      }

      if ($s1 eq "TER")
      {
         ;
      }
      elsif (($s5 >= 1) && ($s5 <= $BIG_NUM))
      {
         if ($DB_WRITE_L) {printf("  write line, s5 num match: $s5, ext: $ext.\n");}
         ## s5 is a residue number. ##
         printf outF ("%-6.6s %4d %s %-4.4s%3.3s %5.4s %11.3f %7.3f %7.3f %5.2f %5.2f \n",
                                              $s1,$s2,"",$s3,$s4,$s5,$s6,$s7,$s8,$s9,$s10);
      }
      else 
      {
         if ($DB_WRITE_L) {printf("  write line, s5, NO num match: $s5, ext: $ext.\n");}
         ## s5 is a substructure designation, s6 will be the residue number. ##
         printf outF ("%-6.6s %4d %s %-4.4s%3.3s %1.1s %4s %10.3f %7.3f %7.3f %5.2f %5.2f \n",
                                              $s1,$s2,"",$s3,$s4,$s5,$s6,$s7,$s8,$s9,$s10,$s11);
      }

      $i++;
   }

   close (outF);
   if ($DB_WRITE) {printf("  returning from WritePDB.\n");}
}