#!/usr/bin/perl -w
use strict ;
#
# opens validate.dat and extracts residue id list
# residue ids have to be unique, multiple chains have to be taken into account (in res id include MOL1)
# write walk_prot.com walks down the protein chain and does sth with each residue
# calls analize_validate.pl for comparing validate.dat modified and original
#
# analize_validate.pl generates sum over each residue and compares old/new
# restores old or keeps new state depending on score
# for now only electron density is being compared in validate.dat
#
#my $density_mode = "gaussian" ; #linear
my $density_mode = "linear" ;
my $file_orig = "validate_starting.dat" ;
open ( FILE_ORIG, "< $file_orig" || die "$! $file_orig" ) ;

my @file_orig = <FILE_ORIG> ;

close ( FILE_ORIG ) ;

sum ( @file_orig ) ;


sub sum
{
    my @file =  @_ ;
    my $line;
    my $nres = 0;
    my @words;
    my @my_seq = ();

    foreach $line ( @file )
    {
	chomp ( $line ) ;
	@words = split(" ", $line) ;
	# words [9-13]
        # in perl one can do array of hashes etc.

	$my_seq[$nres]{"atom_num"} = $words[9] ;     # atom number
	$my_seq[$nres]{"atom_name"} = $words[10] ;   # CA, N, CB
	$my_seq[$nres]{"res_name"} = $words[11] ;    # MET, GLY,
	$my_seq[$nres]{"res_num"} = $words[12] ;     # 1, 2, 
	$my_seq[$nres]{"seg_name"} = $words[13] ;    # MOL1,
	$nres++ ;

    }


# residue is unabigously defined by start/end atom number
# this part determines residue atoms and saves them by atom number

    my @residues = ();
    $residues[0]{"res_name"} = "XXXXXXXX" ;
    $residues[0]{"res_num"} = "XXXXXXXX" ;
    $residues[0]{"seg_name"} = "XXXXXXXX" ;
    my $numb = 0;        # residue count
    my $count = 0;

    for ($count = 0; $count < @my_seq ; $count++ )
    {
	if ( $my_seq[$count]{"res_name"} ne $residues[$numb]{"res_name"}
 	     or $my_seq[$count]{"res_num"} ne $residues[$numb]{"res_num"}
	     or $my_seq[$count]{"seg_name"} ne $residues[$numb]{"seg_name"} )
	{
	    $numb++ ;
	    
	    $residues[$numb]{"res_name"} = $my_seq[$count]{"res_name"} ;
	    $residues[$numb]{"res_num"} = $my_seq[$count]{"res_num"} ;
	    $residues[$numb]{"seg_name"} = $my_seq[$count]{"seg_name"} ;
	    $residues[$numb]{"atom_start"} = $my_seq[$count]{"atom_num"} ;    # starting atom number of residue
	    
#	    print $numb . "\t" . $residues[$numb]{"res_name"} . "\t" . $residues[$numb]{"atom_start"} . "\n";
	}
    }
    shift @residues;	# remove first elem, becase we start writing at $numb = 1

	

# recognize which residue we have and what to do with it
# do stuf and center next
# on si mora zdej zapomnit da sem ze sel skozi kasno ak? ne ni treba
# ampak atom_start bomo rabli za enolicno dolocit preostanek - read res at atom num

# define functions; fit-_side etc.
# make array with all the commands, push with newline

    my @commands = () ;

    my $subroutine = "subroutine char FILE_VALID\n";
    my $SET_WEIGHTS = "set weights select all end = 1.0\n";
    my $FIT_SIDE  = "< ?MAIN_CMDS_NEW:auto_stuff_side.cmds WORK_SEGM FIT\n" ;
    my $RIG_SIDE  = "< ?MAIN_CMDS_NEW:auto_stuff_rig_side.cmds WORK_SEGM\n" ;
    my $FLIP_SIDE = "< ?MAIN_CMDS_NEW:auto_stuff_side.cmds WORK_SEGM FLIP\n";
    my $min = "minimize key\n";    # at some poing smart minimize = minimize until better
    my $min_5 = "minimize key\nminimize key\nminimize key\nminimize key\nminimize key\n";
    my $compare = "SHELL analize_validation_results.pl\n< val.com\n";
    my $center_anal = "anal atom init
ener all off density on
ener density ( scale -1.0 $density_mode )
ener sele active .a .n atom name H* end
anal sele active .a .n atom name H* end energy resid 1
open unit 65 file FILE_VALID over_write
show unit 65 anal sele active .a .n atom name H* C N O CA CB .a weigh 0.01 100. end
clo uni 65
ener density ( scale 10.0 polynom 2 )
ener chemistry on dens on\n";
    

    push @commands, $subroutine;
    push @commands, $SET_WEIGHTS;
#    push @commands, "set level 3\n";
    push @commands, "\n";

    push @commands, "!replaces SHELL create_anal_center.pl
set weights select all end = 1.0
anal atom init
ener all off density on
ener density ( scale -1.0 $density_mode )
ener sele active .a .n atom name H* end
anal sele active .a .n atom name H* end energy resid 1
open unit 65 file validate_starting.dat over_write
show unit 65 anal sele active .a .n atom name H* C N O CA CB .a weigh 0.01 100. end
clo uni 65
ener density ( scale 10.0 polynom 2 )
ener chemistry on dens on

! MAIN part
!key active_starting sele active end
!key active_drops sele active end
!key active sele by res atom numb  end

set coor save now
set coor save off

set vari SAVE_STATUS char = SEL_ACTIVE
set vari SEL_ACTIVE global = XXXX
!key active_drops sele active end\n"; 


    for ($count = 0; $count < @residues; $count++)
    {

	if ( $residues[$count]{"res_name"} eq "GLY" or $residues[$count]{"res_name"} eq "ALA" or $residues[$count]{"res_name"} eq "PRO" )
	{
	    push @commands, "!$count residue is $residues[$count]{'res_name'}, first atom numb is $residues[$count]{'atom_start'}\n";
	    push @commands, "key active select by resid atom numb $residues[$count]{'atom_start'} end\n";
#	    do nothing
	    push @commands, "\n";
	}

	if ( $residues[$count]{"res_name"} eq "MET" or $residues[$count]{"res_name"} eq "CYS" or $residues[$count]{"res_name"} eq "SER" )
	{
	    push @commands, "!$count residue is $residues[$count]{'res_name'}, first atom numb is $residues[$count]{'atom_start'}\n";
	    push @commands, "key active select by resid atom numb $residues[$count]{'atom_start'} end\n";

	    push @commands, $FIT_SIDE;
	    push @commands, $min_5;
	    push @commands, $SET_WEIGHTS;
	    push @commands, $center_anal;
	    push @commands, $compare;
	    push @commands, "\n";
	}
	
	if ( $residues[$count]{"res_name"} eq "ARG" or $residues[$count]{"res_name"} eq "ASP"
	     or $residues[$count]{"res_name"} eq "GLU" or $residues[$count]{"res_name"} eq "LYS" )
	{
	    push @commands, "!$count residue is $residues[$count]{'res_name'}, first atom numb is $residues[$count]{'atom_start'}\n";
	    push @commands, "key active select by resid atom numb $residues[$count]{'atom_start'} end\n";

	    push @commands, $FIT_SIDE;
	    push @commands, $min_5;
	    push @commands, $RIG_SIDE;
	    push @commands, $min_5;
	    push @commands, $SET_WEIGHTS;
	    push @commands, $center_anal;
	    push @commands, $compare;
	    push @commands, "\n";
	}

	if ( $residues[$count]{"res_name"} eq "TYR" or $residues[$count]{"res_name"} eq "PHE"
	     or $residues[$count]{"res_name"} eq "TRP" )
	{
	    push @commands, "!$count residue is $residues[$count]{'res_name'}, first atom numb is $residues[$count]{'atom_start'}\n";
	    push @commands, "key active select by resid atom numb $residues[$count]{'atom_start'} end\n";

	    push @commands, $FIT_SIDE;
	    push @commands, $min_5;
	    push @commands, $RIG_SIDE;
	    push @commands, $min_5;
	    push @commands, $SET_WEIGHTS;
	    push @commands, $center_anal;
	    push @commands, $compare;
	    push @commands, "\n";
	}

	if ( $residues[$count]{"res_name"} eq "VAL" or $residues[$count]{"res_name"} eq "ILE"
	     or $residues[$count]{"res_name"} eq "LEU" or $residues[$count]{"res_name"} eq "THR" )
	{
	    push @commands, "!$count residue is $residues[$count]{'res_name'}, first atom numb is $residues[$count]{'atom_start'}\n";
	    push @commands, "key active select by resid atom numb $residues[$count]{'atom_start'} end\n";

 	    push @commands, $FIT_SIDE;
# 	    push @commands, $min_5;
# 	    push @commands, $RIG_SIDE;
# 	    push @commands, $RIG_SIDE;
# 	    push @commands, $min_5;
# 	    push @commands, $min_5;
# 
 	    push @commands, $SET_WEIGHTS;
	    push @commands, $center_anal;
	    push @commands, $compare;
	    push @commands, "\n";

#	    FLIP_SIDE, min, RIG_SIDE, min, FIT_SIDE, min
#	    little later aligator
	    push @commands, "!$count residue is $residues[$count]{'res_name'}, first atom numb is $residues[$count]{'atom_start'}\n";
	    push @commands, "!FLIP_SIDE\n";
	    push @commands, "key active select by resid atom numb $residues[$count]{'atom_start'} end\n";

	    push @commands, $FLIP_SIDE;  #flip_side
	    push @commands, $min_5;
	    push @commands, $RIG_SIDE;
	    push @commands, $RIG_SIDE;
	    push @commands, $min_5;
	    push @commands, $min_5;
	    push @commands, $FIT_SIDE;   #fit_side
	    push @commands, $min_5;
	    push @commands, $min_5;

	    push @commands, $SET_WEIGHTS;
	    push @commands, $center_anal;
	    push @commands, $compare;
	    push @commands, "\n";

	}
	print "residue name $count <", $residues[$count]{"res_name"}, ">\n" ;
	if ( $residues[$count]{"res_name"} eq "ASN" or $residues[$count]{"res_name"} eq "GLN"
	     or $residues[$count]{"res_name"} eq "HIS" )
	{
	    push @commands, "!$count residue is $residues[$count]{'res_name'}, first atom numb is $residues[$count]{'atom_start'}\n";
	    push @commands, "key active select by resid atom numb $residues[$count]{'atom_start'} end\n";

	    push @commands, $FIT_SIDE;
	    push @commands, $min_5;
	    push @commands, $RIG_SIDE;
	    push @commands, $min_5;
	    push @commands, $SET_WEIGHTS;
	    push @commands, $center_anal;
	    push @commands, $compare;
	    push @commands, "\n";
	}
    }
    push @commands, "<?MAIN_CMDS:re_image.cmds WORK_SEGM\n";
    push @commands, "return\n";
 

    # writes commands to file
    my $file_com = "walk.com" ;
    open ( FILE_COM, "> $file_com" || die "$! $file_com" ) ;

    foreach (@commands)
    {
	print FILE_COM $_ ;
    }
    close ( FILE_COM ) ;
}