#!/usr/bin/env perl
#*******************************************************************************
# Copyright (C) 2001-2004 by Novartis Institutes of Biomedical Research
# 
# This program is free software; you can redistribute it and/or modify it
# under the terms of the Novartis Open Source License as published by
# Novartis AG; either version 1, or (at your option) any later version of
# the license.
# 
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# Novartis Open Source License for more details.
# 
#*******************************************************************************

use POSIX;
use DBI;
use Getopt::Long;
use Pod::Usage;
use Term::ANSIColor;
use Term::ReadLine;
use IO::File;
use File::Basename;
use strict;

my $PRG_PATH=`pwd`;
chomp($PRG_PATH);

my $BASEPATH = getProgPath();

my $HOME  = $ENV{HOME};
my $BIN   = "$BASEPATH/../src";
my $SHARE = "$BASEPATH";

my $user = `whoami`;
chomp($user);

my %mysqlsetup = (
									host => "",
									user => $user,
									passwd => "",
									ask    => 0
								 );

my %expsetup = (
								group => "",
								full => 0,
								path  => $PRG_PATH,
								format => "own",
								probe => 1.4,
								radius => 3.8,
								rotstep => 2.0,
								resolution => 13,
								hyphen => 1,
								script => "$SHARE/run.sh.in"
							 );

# READLINE
my $term = new Term::ReadLine 'Surface Comp';
my $OUT = $term->OUT;
$term->ornaments(0);

#********************************************************************************
# MAIN
#********************************************************************************

my $help;
my $man;

GetOptions('help' => \$help,
					 'man|m' => \$man,
					 'host|h=s' => \$mysqlsetup{host},
					 'user|u=s' => \$mysqlsetup{user},
					 'passwd|p' => \$mysqlsetup{ask},
					 'group|g=s' => \$expsetup{group},
					 'full|f' => \$expsetup{full},
					 'path=s' => \$expsetup{path},
           'format=s' => \$expsetup{format},
					 'probe=f' => \$expsetup{probe},
					 'radius|r=f' => \$expsetup{radius},
					 'rotstep=f' => \$expsetup{rotstep},
					 'resolution=f' => \$expsetup{resolution},
					 'hyphen=i' => \$expsetup{hyphen},
					 'script=s' => \$expsetup{script});


GetOptions() or pod2usage(2);

pod2usage(1) if $help;
pod2usage(-exitstatus => 0, -verbose => 2) if $man;

#-------------------------------------------------------------------------------
# INITIALIZATION

#DIRECTORIES

$expsetup{path} = $PRG_PATH.'/'.$expsetup{path} if ($expsetup{path} !~ /^\//);

my %dirs = checkPath(%expsetup)
 or
	die ("The experiment data path ".
			 "$expsetup{path} does not have then sufficient directory structure!");

my $db = connect2Mysql(%mysqlsetup);

my @surfaces = collectSurfaces(format => $expsetup{format},
															 probe => $expsetup{probe},
															 %dirs);

createDatabase(db => $db, group => $expsetup{group});

foreach my $surface (@surfaces)
{
	readCenters(surface => $surface,
							path => $dirs{center});

	addSurface(surface => $surface,
						 format  => $expsetup{format},
						 db      => $db);
}


my @experiments = createExperiments(db => $db,
																		surfaces => \@surfaces,
																		%expsetup);

createScript(experiments => join(" ", @experiments),
						 %expsetup,
						 %mysqlsetup);

1;

#*******************************************************************************
sub getProgPath
#*******************************************************************************
{
	return dirname($0);
}

#*******************************************************************************
sub checkPath
#*******************************************************************************
{
	my %params = @_;


	my %dirnames = ( surf   => $params{path}."/".$params{format},
									 data   => $params{path}."/data",
									 pdb    => $params{path}."/pdb",
									 center => $params{path}."/center");

	foreach my $dir (values(%dirnames))
	{
		return 0 if (! -d $dir);
	}

	return %dirnames;
}

#*******************************************************************************
sub connect2Mysql
#*******************************************************************************
{
	my %params = @_;

	if ($params{ask})
	{
		print($OUT "Give a password for ".
					"$params{db} on $params{host} as user $params{user}:\t");
		system("stty -echo"); # Do not echo the password at the prompt!!!
		$params{passwd} = <STDIN>;
		system("stty echo");
		print($OUT "\n");
	}

	my $data_source="DBI:mysql:host=" . $params{host};
	my $dbh = DBI->connect($data_source,$params{user},$params{passwd})
		or die("Could not connect to database $params{db} on host $params{host}");

	return $dbh;
}

#*******************************************************************************
sub createDatabase
#*******************************************************************************
{
	my %params = @_;

	my $create = $params{db}->prepare("CREATE DATABASE $params{group};");
	$create->execute() or die($create->errstr);

	$params{db}->do("use $params{group}");

	open(DBDEFS, "$SHARE/dbdefinition.sql");
	my @commands = split(/\n\n/,join("", <DBDEFS>));
	close(DBDEFS);
	
	foreach my $command (@commands)
	{
		$params{db}->do($command);
	}
}


#*******************************************************************************
sub collectSurfaces
#*******************************************************************************
{
	my %params = @_;
	my @surfaces;

	opendir(DIR, $params{surf}) || die "can't opendir $params{surf}: $!";
	my @files = grep { /\.$params{format}$/ && -f "$params{surf}/$_" } readdir(DIR);
	closedir DIR;

	@files = sort(@files);

	foreach my $file (@files)
	{
		my $name = $file;
		$name =~ s/\.$params{format}$//;

		my $pdb  = "$params{pdb}/$name.pdb";
		my $data = "$params{data}/$name.dat";
		
		(-f $pdb) or
			die("No structure file for surface $name");
		(-f $data) or
			die("No data file for surface $name");

		$file = "$params{surf}/$file";

		#retrieve necessary surface information
		my $surfinfo = `$BIN/surfinfo -i $file -f $params{format} points area`;
		chomp($surfinfo);
		my ($points, $area) = split(/\n/,$surfinfo);
		$area=$points/3 if ($area==0);

		push(@surfaces, { file => $file,
											name => $name,
											pdb => $pdb,
											data => $data,
											density => $points/$area,
											probe => $params{probe},
											points => $points,
										});
	}

	return @surfaces;
}

#*******************************************************************************
sub readCenters
#*******************************************************************************
{
	my %params = @_;
	my $surface = $params{surface};
	my $centerfile = "$params{path}/$surface->{name}.txt";

	my %centers;

	if (-f $centerfile)
	{
		open(CENTERFILE, "<$centerfile");
		while (<CENTERFILE>)
		{
			chomp;
			s/^\s+//;
			my($x, $y, $z, $r, $spec) = split(/\s+/);
			$centers{$spec} = { x => $x, 'y' => $y, z => $z, r => $r};
		}
		close(CENTERFILE);
	}

	$surface->{centers} = { %centers };
}


#*******************************************************************************
sub addSurface
#*******************************************************************************
{
	my %params = @_;
	my $surface = $params{surface};
	my $format  = $params{format};
	my $db      = $params{db};

	my $surfcode = $db->quote(convertSurface(file   => $surface->{file},
																					 format => $format));


	my @vals = ( 0,
							 "'$surface->{name}'",
							 "'$surface->{pdb}'",
							 $surface->{probe},
							 $surface->{density},
							 $surface->{points},
							 $surfcode);
	my $val = join(", ", @vals);
	my $insert = $db->prepare("INSERT INTO surfaces VALUES ($val);");

	$insert->execute() or die($insert->errstr);
	$surface->{id} = $insert->{mysql_insertid};

	my $center_insert = 
		$db->prepare("INSERT INTO centers VALUES (0, ?, ?, ?, ?, ?);");
	my $cs_insert = 
		$db->prepare("INSERT INTO surface_centers VALUES ($surface->{id}, ?);");
	
	foreach my $spec (keys(%{$surface->{centers}}))
	{
		my $center = $surface->{centers}->{$spec};
		$center_insert->execute($center->{x},
														$center->{'y'},
														$center->{z},
														$center->{r},
														$spec);

		$center->{id} = $center_insert->{mysql_insertid};
		$cs_insert->execute($center->{id});
	}
}

#*******************************************************************************
sub convertSurface
#*******************************************************************************
{
	my %params = @_;
	my $file   = $params{file};
	my $format = $params{format};

	my $code = "";
	
	system("$BIN/surfaceconvert -i $file -o temp.bin $format bin");
	
	open(SURFACE, "<temp.bin");
	binmode(SURFACE);
	
	my $buffer;
	while (read(SURFACE, $buffer, 1024))
	{
		$code .= $buffer;
	}

	close(SURFACE);
	unlink("temp.bin");

	return $code;
}

#*******************************************************************************
sub createExperiments
#*******************************************************************************
{
	my %params = @_;
	my $db       = $params{db};
	my $surfaces = $params{surfaces};
	my $full     = $params{full};
	my $hyphen   = $params{hyphen};

	my @experiments;

	my $expsql = $db->prepare("INSERT INTO experiments VALUES(".
													 " 0, ?, $params{radius}, 1, ?, ?,".
													 " $params{resolution}, $params{rotstep},".
													 " ?, ?, 0, '', ?, ?);");

	#LOOP for surface A
	for(my $i=0;$i<=$#{$surfaces};++$i)
	{
		my $A = $$surfaces[$i];

		#LOOP for surface B
		for(my $j=$i+1;$j<=$#{$surfaces};++$j)
		{
			my $B = $$surfaces[$j];
			
			#LOOP for common sites
			foreach my $spec (keys(%{$A->{centers}}))
			{
				if (exists($B->{centers}->{$spec}))
				{
					$expsql->execute(createName(a => $A->{name}, b => $B->{name},
																			c => $spec, hyphen => $hyphen),
													 $A->{id},
													 $B->{id},
													 $A->{data},
													 $B->{data},
													 $A->{centers}->{$spec}->{id},
													 $B->{centers}->{$spec}->{id});
				
					push(@experiments, $expsql->{mysql_insertid});
					
				}
			}
		
			#add full comparison if necessary
			if ($full)
			{
				$expsql->execute(createName(a => $A->{name}, b => $B->{name},
																		hyphen => $hyphen),
												 $A->{id},
												 $B->{id},
												 $A->{data},
												 $B->{data},
												 0, 0);
			
				push(@experiments, $expsql->{mysql_insertid});
			}
		}
	}

	return @experiments;
}

#*******************************************************************************
sub createName
#*******************************************************************************
{
	my %params = @_;
	my $hyphen = $params{hyphen};
	
	my @a = split(/\-/, $params{a});
	my @b = split(/\-/, $params{b});

	if (--$hyphen>=0)
	{
		pop(@a) while ($#a>($hyphen));
		pop(@b) while ($#b>($hyphen));
	}
	
	my @name = (@a, @b);
	push(@name, $params{c})if (exists($params{c}));

	return join("-", @name);
}

#*******************************************************************************
sub createScript
#*******************************************************************************
{
	my %params = @_;
	
	my %replaces = ( DBHOST      => $params{host},
									 DB          => $params{group},
									 EXPERIMENTS => $params{experiments});

	my $scripttempl = $params{script};
	my @scriptparts = split(/\//, $scripttempl);
	my $scriptname  = pop(@scriptparts);
	$scriptname =~ s/\.in$//;

	open(INFILE,  "<$scripttempl");
	open(OUTFILE, ">$scriptname");

	while (<INFILE>)
	{
		foreach my $replace (keys(%replaces))
		{
			s/\@$replace\@/$replaces{$replace}/;
		}
		print OUTFILE;
	}

	close(INFILE);
	close(OUTFILE);

	chmod(0755, $scriptname);
}

__END__

=head1 NAME

preparesurfcomp

=head1 SYNOPSIS

A script that prepares SURFCOMP experiments

B<preparesurfcomp [ --help --man|m ]>
                  B<--host|h> I<DATABASE_HOST> B<--user|u> I<DB_USER> B<--passwd|p>
                  B<--group|g> I<GROUPNAME> B<--path> I<DATAPATH> B<--full|f>
                  B<--format> I<FILE_FORMAT>
                  B<--probe> I<PROBE_RADIUS> B<--radius|r> I<PATCH_RADIUS>
                  B<--rotstep> I<ROTATION_STEP> B<--resolution> <HSI_RESOLUTION>
                  B<--hyphen> I<#NAME_PARTS> B<--script|s> I<RUNSCRIPT>

=head1 OPTIONS

=over 4

=item --help

print out this help page.

=item --man

print out a detailed help page as manpage

=item --host|h [default=I<current machine>]

The MYSQL database host where the experimental database should be created and
stored. You can either use short names or full qualified names.

=item --user|u [default=I<current user>]

The name of the database user that should be used for the DB operations. If the
user has a password you have to give the option B<-p> or B<--passwd> as well.

=item --passwd|p

If this flag is specified the script prompts for the MYSQL password of the
given user

=item --group|g

The name of the group of experiments that will be prepared by the script. The
groupname will be used as the name of the experimental database.

=item --path [default=I<.>]

The path of the input data for the experiments. (see also
L</"INPUT DATA STRUCTURE">)

=item --full|f

If this flag is specified a full comparison experiment (not regarding any
surface sites) is created for each pair of surfaces.

=item --format [default=I<own>]

The file format of the surface files. Possible values are off, I<own>, msms and
sld. Please note that the specification of the format determines the name of
the directory where the surface files are expected.

=item --probe [default=I<1.4 Angs.>]

The size of the probe sphere that was used to create the molecular
surfaces. This parameter is for information purposes only and can be left as is
or set to any value if no solvent accessible surfaces have been used.

=item --radius|r [default=I<3.8 Angs.>]

The radius of the surface patches. Values between 3.5 and 6 Angs. are
appropriate.

=item --rotstep [default=I<2°>]

The size of an angular rotation step in the harmonic shape image
filtering. This parameter need not be specified in almost any cases.

=item --resolution [default=I<13>]

Lateral resolution of the harmonic shape images. Any value between 11 and 15
will be fine.

=item --hyphen [default=I<1>]

Determines how many parts of the surface names will be used for the experiment
names (see below).

=item --script|s [default=I<(datadir)/run.sh.in>]

The script template for the generation of then surfcomp run script.

=back

=head1 DESCRIPTION

The B<preparesurfcomp> script builds the experimental database for the
comparison of a set of surfaces and surface sites which are stored in a given
input data directory. It can also generate a shell script that can invoke
B<surfcomp> for all experiments. It therefore simplifies to a large extend the
preparation of B<surfcomp> computations.

Although a lot of parameters can be specified, those that have only a long name
can usually be omitted. And most of the other parameters are have quite
reasonable defaults.

=head2 INPUT DATA STRUCTURE

The script expects surface files, surface point data files, structure
files and maybe center descriptions. All these files have to be stored in the
correct directory structure:

   /
   +-center           (the surface site sphere descriptions)
   +-data             (the surface point data files)
   +-pdb              (the structure file for each surface)
   +-own|msms|sld|off (the surface files in the correct format)

=head2 SCRIPT TEMPLATE

The final step of the preparation is the generation of a run script that
automatically invokes B<surfcomp> for all generated experiments. This is
usually a shell script (Bourne or C-shell) but can as well be written in any
other scripting language. B<preparesurfcomp> uses a template file to generate
this script. Such a template file contains all the commands and parameter
definitions that are necessary for the execution of the B<surfcomp> program.

In this template three placeholders will be replaced by the appropriate values
for the generated group of experiments:

=over 4

=item @DBHOST@

will be replaced by the name of the database host where the experiment database
is stored.

=item @DB@

the name of the experiment database (actually the group name).

=item @EXPERIMENTS@

the list of experiment IDs that have been created in the current group.

=back

A special task of the run script is the definition of the various filter
thresholds and data columns. Please referre to L<surfcomp(1)> for the exact
specfications of those parameters.

The default installation includes a basic Bourne- and C-shell script which
invoke surcomp on the local machine with a basic set of filter heuristics I use
usually for the comparison of small molecular surfaces. Any other script
template can be defined and passed to the program by the B<--script|s>
parameter.
