#!/usr/bin/perl
#Author: Le Quy Vang
#Purpose: +This script will extract atom names and their coordinates 
#         and save them into files given by users. It then calculate
#         longest distance from the atom input by user to the rest of
#         atom in the input data file. 
#        +Moreover, this script can use directly PDB file as input data 
#         file. It also can offer user opportunity to input other atom 
#         or data file for calculation. 
#         
#Date: November, 2006
use strict;

print "This program find maximum distance between a given atom and \n";
print "another atom among the remaining ones \n";

my $OtherFile;
do {
    
    fileopen(); # this is a subroutin call for opening file;
    
    #Get interested content
    my $AtomNameList='';
    my $AtomCoordinates='';
    foreach (<INPUTFILE>)
     {
	    my $line1 =substr($_,13,13); 
	    my $line2 =substr($_,31,25);
	   $AtomNameList.= $line1."\n"; 
	   $AtomCoordinates.= $line2."\n"; 
      }
    close INPUTFILE;
    
    #Write to out file
    print "File was successfully read \n";
    print "Do you want to get atom list and coordinate files ? Y/N \n";
    
    #Get answer from user about exporting file for further use.
    my $answer = <STDIN>;
    if (uc($answer)=~ m/Y/)
      {
        SaveExtracted($AtomNameList,$AtomCoordinates);
      }
    
    # Call subroutine to get arrays of Atom Name and Atom Coordinates respectively   
	my @arrAtomName= BreakDown2Arr($AtomNameList,"\n");
	my @arrAtomCoordinates =BreakDown2Arr($AtomCoordinates,"\n");
	   
   my $ContinueOtherAtom; #This variable to catch user answer for continuing on other atom
  do
     {
      #Get input for calculation from user
     print "Which atom do you want to calculate maximum distance?\n";
     print "Please enter as following format: Residuename:ResidueNo:Atomname\n";
     my $intAtomRow=-1;
     my $atom='';
    do 
     {   
           if ($atom ne '')
            {
                print "Atom could not be found, please re-enter\n";
            }

              $atom=<STDIN>;
              chomp($atom);
              $atom=uc($atom);
	      if ($atom !~/[A-Z]{3}:[1-9]{1,4}:[A-Z]{1,2}[1-9]{0,2}/)
              {
                   print "Invalid input format\n";
                   print "Input should be like ASP:113:OD1 \n";    
              }
	      $intAtomRow=FindAtom($atom,@arrAtomName); 
     } while ($intAtomRow<0);
     
      #Calculate distance array and pin point the maximum value:
      my @arrDist=();
      my $DistMax=0;
      for (my $i = 0; $i < @arrAtomCoordinates; $i ++)
          {
             my $dist=sqrt(
                           (($arrAtomCoordinates[$i]->[0])-($arrAtomCoordinates[$intAtomRow]->[0]))**2 +
                           (($arrAtomCoordinates[$i]->[1])-($arrAtomCoordinates[$intAtomRow]->[1]))**2 +
                           (($arrAtomCoordinates[$i]->[2])-($arrAtomCoordinates[$intAtomRow]->[2]))**2
                          );
	     $arrDist[$i]=$dist;
              if ($dist > $DistMax)
               {
                $DistMax = $dist;
                
               }
          }
      # There might be multiple atom have longest distance! Try to catch all of them here.

       my @MaxAtomArr=();
       for (my $i=0;$i<@arrDist;$i++)
       {
	if ($arrDist[$i] == $DistMax)
	 {
         push (@MaxAtomArr,$i);
    	 }
       }
       #print out result. 
       print "\n MAX DISTANCE =", $DistMax, "\n";
       foreach my $i (@MaxAtomArr)
       {
         my $strAtomMax= ($arrAtomName[$i]->[1].$arrAtomName[$i]->[2].":". $arrAtomName[$i]->[0]);
           
           print "The maximum distance atom is at row:", $i+1, "\n";
           print "Atom name is ", $strAtomMax , "\n";
           print "Position of $strAtomMax is at  ", $arrAtomCoordinates[$i]->[0], "\t",
                                                      $arrAtomCoordinates[$i]->[1], "\t",
                                                      $arrAtomCoordinates[$i]->[2], "\n";
       }
    


         print "Do you want to calculate other atom? Y/N \n";
         $ContinueOtherAtom =<>;
          
   } while ($ContinueOtherAtom !~ m/n|N/);
     print "Do you want to work with other file? Y/N, \n";
       $OtherFile=<>;
} while ($OtherFile !~ m/n|N/);
   print "Goodbye! \n";
#
# Main program ends!   
#
#
# This section start subroutines' code:   
#   
      #----------------------------------------------------------------------------------------
  #Sub fileopen
  #This sub is to looping for opening input file
    sub fileopen
    {
     #Get file to process
      my $FileIn='';
      
    do {   #do untill get valid file to open.
	
	if ($FileIn ne '')
	{      
	       print "File could not be found! \n";
	       print "Please re-enter filename: \n";
	}
	print "Please specify input filename, file MUST be in PDB format: \n";
	
	$FileIn=<STDIN>;
	 chomp($FileIn);
    }while (!open  INPUTFILE, $FileIn);
    return
    }
     #----------------------------------------------------------------------------------------   
    #Subroutin Savefile
    sub SaveExtracted
    {
           my ($subNameList, $subCoordinates)=@_;      
      print "Please give filename for Atom List: \n";
           my $fileAtomNameList=<>;
	   open ANAMEOUT, ">$fileAtomNameList" or die 'Can\'t open the file for writing';
	   print ANAMEOUT $subNameList;
	   close ANAMEOUT;
          print "Please give filename for Coordinate List: \n";
          my $fileAtomCoordinates=<>;
	   open NXYZOUT, ">$fileAtomCoordinates" or die 'Can\'t open the file for writing';
	   print NXYZOUT $subCoordinates;
	   close NXYZOUT;
          # print "Data files have been saved \n";
    }
    #----------------------------------------------------------------------------------------
     #Subroutine to break input content into an array.
    sub BreakDown2Arr 
    {
        my ($substrB4Break, $rowbreak)=@_;
	my @arrBreak= split(/$rowbreak/, $substrB4Break);
	my @arrAfter=();
	my $rowCount=0;
     foreach (@arrBreak)
     {
	    my @vars = split;
	    for (my $i = 0; $i < @vars; $i ++)
	    {
	     $arrAfter[$rowCount]->[$i] = $vars[$i];
	    }
	    $rowCount++;
      }    
       return @arrAfter;
    }
    #----------------------------------------------------------------------------------------
    sub FindAtom
    {
     my ($toFind, @arrIn)=@_;
     my ($resid, $resno, $at)='';
     my $RowNum=-1;
     
              ($resid,$resno,$at)=split(/:/,$toFind);
              for (my $i=0; $i<(@arrIn); $i++) # Search for input atom in data array. 
              {
                   if (($arrIn[$i]->[2] == $resno)&&($arrIn[$i]->[1] eq $resid)&&($arrIn[$i]->[0] eq $at))
                   {
                      $RowNum=$i;
                      
                   }
              }
      return $RowNum;
    }
#
# subroutines' code ends  
#