#!/usr/bin/perl -w
use strict;
use warnings;

####Genetic Drift###########
print "#######Testing Genetic drift effects with differing sample sizes######\n\n";
print "Set the parameters below:\nPopulation Size (between 100 and 100000): \n";
my($pop_size, $temperature, $oxygen, $food, $mutation_rate, $number_generations, $bottleneck, $initial_organism, $student_number);
$pop_size = <STDIN>;
print "Number of generations (between 20 and a 1000): \n";
if(($pop_size <1) || ($pop_size >1000))
  {
    print "What if you try to use a sensible size!!\n";
    exit;
  }
$number_generations = <STDIN>;
chomp $number_generations;
if($number_generations < 2)
  {
    print "You don't seem up for multiplying do you :D\nTry again, don't be lazy :D\n\n";
    exit;
  }
print "\nTemperature of your environment (1 to 9) (1: cold; 9: hot):";
$temperature = <STDIN>;
chomp $temperature;
if(($temperature < 1) || ($temperature > 9))
  {
    print "\nYea yea, very funny but you have to follow the rules!!!\nTry again ;)\n\n";
    exit;
  }
print "\nOxygen in the air (1 to 9) (1: little oxygen; 9: very oxygenized):";
$oxygen = <STDIN>;
chomp $oxygen;
if(($oxygen < 1) || ($oxygen > 9))
  {
    print "\nYea yea, very funny but you have to follow the rules!!!\nTry again ;)\n\n";
    exit;
  }
print "\nAmount of available food (1 to 9) (1: very scarse; 9: abundant):";
$food = <STDIN>;
chomp $food;
if(($food < 1) || ($food > 9))
  {
    print "\nYea yea, very funny but you have to follow the rules!!!\nTry again ;)\n\n";
    exit;
  }
print "Mutation rate (between 0.1 and 0.9): \n";
$mutation_rate = <STDIN>;
chomp $mutation_rate;
print "Bottleneck to population size (between 0.1 and 0.9): \n";
$bottleneck = <STDIN>;
chomp $bottleneck;
if(($bottleneck < 0) || ($bottleneck > 1))
{
  print "Bad choice!! Bottleneck should be a number ranging between 0 and 1\nTry again\n\n";
  exit;
}

print "\nInitial organism:";
$initial_organism = <STDIN>;
chomp $initial_organism;
if(length($initial_organism) == 2)
  {
    my(@a_org);
    $a_org[0] = substr($initial_organism, 0, 1);
    $a_org[1] = substr($initial_organism, 1, 1);
    if(($a_org[0] <= 0) && ($a_org[1] <= 0))
      {
	print "Bad choice!! Organism should be composed of two digits, each ranging between 1 and 9\nTry again\n\n";
	exit;
      }
  }else
  {
    print "Bad choice!! Organism should be composed of two digits, each ranging between 1 and 9\nTry again\n\n";
    exit;
  }
my($prediction);
print "\nWhat is your predicted organism? ";
$prediction = <STDIN>;
chomp $prediction;
if(length($prediction) == 2)
  {
    my(@a_or);
    $a_or[0] = substr($prediction, 0, 1);
    $a_or[1] = substr($prediction, 1, 1);
    if(($a_or[0] <= 0) && ($a_or[1] <= 0))
      {
	print "Bad choice!! Organism should be composed of two digits, each ranging between 1 and 9\nTry again\n\n";
	exit;
      }
  }else
  {
    print "Bad choice!! Organism should be composed of two digits, each ranging between 1 and 9\nTry again\n\n";
    exit;
  }
print "\nIntroduce your student number: ";
$student_number = <STDIN>;
chomp $student_number;
my($out) = "Results_GeneticDrift" . "_$student_number";
print "\n*****Your evolved population will be sent to $out*****\n";
open(O, ">$out");
srand(time|$$);
#my(@frequencies) = (0.05, 0.1, 0.15, 0.2, 0.25);
my($i, $j) = (0,0);
my(@population_values, @population);
#while($i < scalar(@frequencies))
#  {
 #   my($num) = $frequencies[$i] * $pop_size;
 #   my($valA) = $i + 1; my($valB) = 9 - $i;
 #   my($k) = 0;
 #   for($k = 0; $k < $num; $k++)
 #     {
#	$population_values[$j] = $valA;
#	$j += 1;
#	$population_values[$j] = $valB;
#	$j += 1;
 #     }
 #   $i += 1;
 # }
while($i < $pop_size)
  {
  #  $population_values[$i] = $initial_organism;
    $population[$i] = $initial_organism;
     $i += 1;
  }
$i = 0;
#my(@population);
#while($i < $pop_size)
#  {
#    my($body_size, $move);
#    $body_size = $sorted_population[int(rand(scalar(@sorted_population)))];
#    $move =  $sorted_population[int(rand(scalar(@sorted_population)))];
#    $population[$i] = $body_size . $move;
#    $i += 1;
#  }
###We are starting with the evolutionary dynamic now with the mutation process
my($generations) = 0;
my(@population_grown);
my(%Frequencies);
while($generations < $number_generations)
  {
    print "\rGeneration: " . ($generations + 1);
    my($number_mutations) = 0.1 * ($pop_size);
    my(@chosen_organisms);
    $i = 0;
    while($i < $number_mutations)
      {
	my($orden, $sense);
	$orden = int(rand(2));
	$sense = int(rand(2));
	if($sense < 2)
	  {
	    $sense = 1;
	  }else
	    {
	      $sense = -1;
	    }
	my($organism) = int(rand(scalar(@population)));
	if($orden < 2)
	  {
	    $orden = 0;
	  }else
	    {
	      $orden = 1;
	    }
	my($first, $second);
	$first = substr($population[$organism], 0, 1);
	$second = substr($population[$organism], 1, 1);
	if($orden == 0)
	  {
	    if(($first > 1) && ($first < 9))
	      {
		$first += $sense;
	      }elsif($first == 1)
		{
		  $first += 1;
		}elsif($first == 9)
		  {
		    $first -= 1;
		  }
	    $population[$organism] = $first . $second;
	  }else
	    {
	      if(($second > 1) && ($second < 9))
		{
		  $second += $sense;
		}elsif($second == 1)
		  {
		    $second += 1;
		  }elsif($second == 9)
		    {
		      $second -= 1;
		    }
	      $population[$organism] = $first . $second;
	    }
	$i += 1;
      }
    #########Replicating the population############
    my($num_individuals) = 0;
    my(@new_population);
     #########Now the genetic drift process will start here############
    my($population_transfer_next_generation) = $pop_size - ($bottleneck*$pop_size);
    my($pop) = 0;
    while($pop < $population_transfer_next_generation)
      {
	$new_population[$pop] = $population[int(rand(scalar(@population) - 1))];
	$pop += 1;
      }
    
    #########Founder Effect############
    while($pop < $pop_size)
      {
	$new_population[$pop] = $population[int(rand(scalar(@new_population) - 1))];
	$pop += 1;
      }
    my($indis) = 0;
    foreach my $indi(@new_population)
      {
	$population[$indis] = $indi;
      }
    foreach my $freq(@population)
      {
	$Frequencies{$freq} = 0;
      }
    foreach my $fre(@population)
      {
	$Frequencies{$fre} += 1;
      }
    foreach my $FR(keys(%Frequencies))  
      {
	$Frequencies{$FR} /= (scalar(@population));
	$Frequencies{$FR} *= 100;
      }
    #########Frequencies in that generationt############
    print "\n";
    print O "Student: $student_number\n------------------------\n\n";
    print O "Predicted Organism: $prediction\n";
    print O "Initial organism: $initial_organism\nPopulation_size: $pop_size\nMutation Rates: $mutation_rate\nBottleneck to population size: $bottleneck\n";
    print O "\nEnvironment Settings\n---------------------\n";
    if($temperature < 3)
      {
	print O "Temperature ($temperature): Cold\n";
      }elsif(($temperature >= 4) && ($temperature < 7))
	{
	  print O "Temperature ($temperature): mild\n";
	}elsif($temperature >=7)
	  {
	    print O "Temperature ($temperature): Hot\n";
	  }
    if($oxygen < 3)
      {
	print O "Oxygen ($oxygen): rare\n";
      }elsif(($oxygen >= 4) && ($oxygen < 7))
	{
	  print O "Oxygen ($oxygen): oxygen based atmosphere\n";
	}elsif($oxygen >=7)
	  {
	    print O "Oxygen ($oxygen): High oxygen pressure\n";
	  }
    if($food < 3)
      {
	print O "Food resources ($food): scarse\n";
      }elsif(($food >= 4) && ($food < 7))
	{
	  print O "Food resources  ($food): medium\n";
	}elsif($food >=7)
	  {
	    print O "food  ($food): optimal environment with excess of resources\n";
	  }
    print O "\nFinal Population Frequencies:\n--------------------------------\n\n";
    print O "Phenotype\t\tFrequency in the Population (%)\n---------\t\t--------------------------\n";
    foreach my $Fr(sort{$Frequencies{$a} <=> $Frequencies{$b}} (keys(%Frequencies)))
      {
	print O "$Fr\t\t\t$Frequencies{$Fr}%\n";
      }

    $generations += 1
  }

exit;