#!/usr/local/bin/perl --  # -*-Perl-*-

# $Header: /nfs/pathsoft/repository/psu/genexec/perl/src/repeat_finder.pl,v 1.9 1999/10/06 15:35:12 kmr Exp $

# purpose: this script will search for repeats in the sequence for an
# purpose: embl entry and generate a feature for each repeat

use strict;

use Embl;
use EMBL::Parser;


sub usage
{
  die <<EOF;
usage $0 -h (for help) 
   or $0 [-b] [-s] [-t] embl_file

The -b flag will look for matches on both strands.
The -s flag will cause single base repeats to be ignored
The -t flag will cause one and two base repeats to be ignored
EOF
}

if (@ARGV == 0 || $ARGV[0] eq "-h") {
  usage ();
}

my $both = 0;
my $single_base_repeats = 1;
my $two_base_repeats = 1;

if ($ARGV[0] eq '-b') {
  shift;
  $both = 1;
}

if ($ARGV[0] eq '-s') {
  shift;
  $single_base_repeats = 0;
}

if ($ARGV[0] eq '-t') {
  shift;
  $two_base_repeats = 0;
}

open IN_FILE, "$ARGV[0]" or die "$0: can't open $ARGV[0]\n";

# read in some embl files from a stream
my $set = Read Embl(\*IN_FILE); # read all of the embl entries from that stream

# get an array of all of the EMBL::Entry objects within $set
my @entries = $set->entries();

if (@entries != 1) {
  die "$0: $_[0] contains more than one entry\n";
}

my $entry = pop @entries;

my $sequence = $entry->SQ ();
  
my $sequence_length = $sequence->length ();

my $sequence_string = $sequence->seq;

sub write_match
{
  my ($whole_match_length, $match) = @_;

  my $start = (pos $sequence_string) - $whole_match_length + 1;
  my $end = $start + $whole_match_length - 1;
  my $match_count = 1.0 * $whole_match_length / (length $match);

  print <<EOF;
FT   repeat_region   $start..$end
FT                   /note="($match)$match_count"
FT                   /label={$match}$match_count
EOF

  if ($both) {
    $match =~ tr/atgc/tacg/;
    $match = reverse $match;
  print <<EOF;
FT   repeat_region   complement($start..$end)
FT                   /note="($match)$match_count"
FT                   /label={$match}$match_count
EOF
  }
}

for (my $i = 0 ; $i < $sequence_length ; ++$i) {
  
}

while ($sequence_string =~ /((.)\2{7,}|(..)\3{3,}|(...)\4{3,}|(.{4,100})\5{2,})/g) {
  my $match;
  if (defined $2) {
    $match = $2;
    if (!$single_base_repeats || !$two_base_repeats) {
      next;
    }
  }
  if (defined $3) {
    $match = $3;
    if (!$two_base_repeats) {
      next;
    }
  }
  $match = $4 if (defined $4);
  $match = $5 if (defined $5);
  write_match (length $1, $match);
}
