#!/usr/bin/perl
#Author: Feng Cheng
use warnings;
use strict;
use SVG;
#use GD;
#use GD::Arrow;
#use GD::SVG;

my $synG = $ARGV[0];       # b2a_colinear_orthologous_unique.r
my $synF = $ARGV[1];       # final synteny fragments
my $pic = $ARGV[1].".svg";

&processSyn($synG,$synG.".re");

my ($widtn,$height) = (1600,1200);
&process($synG.".re",$synF,$pic,$widtn,$height);
system("rm $synG.re");

sub process {
  my ($synG,$synF,$pic,$widtn,$height) = @_;

  my ($xw,$yh,$x0,$y0,$im,$x1,$y1,$x2,$y2,@color,$maxx,$maxy);
  &setCanvas($widtn,$height,\$xw,\$yh,\$x0,\$y0,\$im);
  &drawAxis($xw,$yh,$x0,$y0,$x1,$y1,$x2,$y2,\$im);

  my ($gposb,$gposa) = ($synG.".b.gene",$synG.".a.gene");
  &getGenePos($synG,$gposb,$gposa);
  my ($bgene2agene) = &readSyn1($synG);
  my (%gene2pos,%index2vchr,@vchr,%index2hchr,@hchr);
  &readpos($gposb,\%gene2pos,\%index2vchr,\@vchr,\$maxy);
  &readpos($gposa,\%gene2pos,\%index2hchr,\@hchr,\$maxx);

  my ($fposb,$fposa) = ($synF.".b.gene",$synF.".a.gene");
  &getGenePos($synF,$fposb,$fposa);
  my ($bf2af) = &readSyn1($synF);
  my (%frag2pos,%f2sites,%f2pos);
  &readposF($fposb,\%f2pos,\%f2sites,\@vchr);
  &readposF($fposa,\%f2pos,\%f2sites,\@hchr);

  &drawXYpos($xw,$yh,$x0,$y0,$x1,$y1,$x2,$y2,\$im,$maxx,$maxy);
  &csep($xw,$yh,$x0,$y0,$x1,$y1,$x2,$y2,\$im,$maxx,$maxy,\@vchr,\@hchr,\%index2vchr,\%index2hchr);

  my $color = "red";
  &drawData($xw,$yh,$x0,$y0,$x1,$y1,$x2,$y2,\$im,$maxx,$maxy,$bgene2agene,\%gene2pos,\@vchr,\@hchr,$color);

  $color = "green";
  &drawDataF($xw,$yh,$x0,$y0,$x1,$y1,$x2,$y2,\$im,$maxx,$maxy,$bf2af,\%f2pos,\%f2sites,$color);

  open(PNG,">$pic");
  print PNG $im->xmlify;
  close(PNG);

  #system("/home/chengf/bin/distributing_svg_4.74/svg2xxx_release/svg2xxx $pic -t png");
  system("rm $gposb $gposa $fposb $fposa");

}

sub getGenePos {
  my ($synG,$gposb,$gposa) = @_;

  system("awk '{print \$1\"\t\"\$2\"\t\"\$3\"\t\"\$4}' $synG > $gposb");
  system("sort -k2,2 -k3,3n $gposb > $gposb.sort");
  system("mv $gposb.sort $gposb");
  system("awk '{print \$5\"\t\"\$6\"\t\"\$7\"\t\"\$8}' $synG > $gposa");
  system("sort -k2,2 -k3,3n $gposa > $gposa.sort");
  system("mv $gposa.sort $gposa");
}

sub readSyn1 {
  my ($in) = @_;
  my %bgene2agene;
  open(FR,$in);
  while(<FR>) {
    #if(!(/Scaffold/) && !(/scaffold_\d\d+/) && !(/NW/)) {       # it depends---
    if(!(/Scaffold/)) {       # it depends---
      my @temp = split(/\t/);
      $bgene2agene{$temp[0]} = $temp[4];
    }
  }
  close(FR);
  return \%bgene2agene;
}

sub readposF {
  my ($in1,$gene2pos,$gene2sites,$vchr) = @_;
  my @tchr;
  open(FR1,$in1);
  my $out1 = $in1.".sort";
  open(FW1,">$out1");
  my @input1 = map{[$_,(split(/\t/))[1],(split(/\t/))[2],(split(/\t/))[3]]} <FR1>;
  my @input2 = sort{$a->[1] cmp $b->[1] || $a->[2] <=> $b->[2] || $a->[3] <=> $b->[3]} @input1;
  my @sorted = map{$_->[0]} @input2;
  for(@sorted) {
    print FW1 "$_";
  }
  close(FR1);
  close(FW1);

  open(FR1,$out1);
  my ($chr,$gene);
  my ($count,$index) = (0,0);
  my @temp;
  while($_=<FR1>) {
    #if(!(/Scaffold/) && !(/^scaffold_\d+/)) {
    if(!(/Scaffold/)) {       # it depends---
      chomp;
      @temp = split(/\t/);
      $gene2pos->{$temp[0]} = ($temp[2]+$temp[3])/2;
      $gene2sites->{$temp[0]} = "$temp[2]\t$temp[3]";
      if($count!=0) {
        if($temp[1] ne $chr) {
          $index += 1;
        }
      }
      else {
        $index += 1;
      }
      $gene2pos->{$temp[0]} += $vchr->[$index-1];
      $count += 1;
      $gene = $temp[0];
      $chr = $temp[1];
    }
  }
  $tchr[$index] = $gene2pos->{$gene} + abs($temp[2]-$temp[3])/2;
  close(FR1);

  system("rm $out1");
}

sub readpos {
  my ($in1,$gene2pos,$index2vchr,$vchr,$maxy) = @_;

  open(FR1,$in1);
  my $out1 = $in1.".sort";
  open(FW1,">$out1");
  my @input1 = map{[$_,(split(/\t/))[1],(split(/\t/))[2],(split(/\t/))[3]]} <FR1>;
  my @input2 = sort{$a->[1] cmp $b->[1] || $a->[2] <=> $b->[2] || $a->[3] <=> $b->[3]} @input1;
  my @sorted = map{$_->[0]} @input2;
  for(@sorted) {
    print FW1 "$_";
  }
  close(FR1);
  close(FW1);

  open(FR1,$out1);
  my ($chr,$gene);
  my ($count,$index) = (0,0);
  while($_=<FR1>) {
    #if(!(/Scaffold/) && !(/scaffold_\d\d+/) && !/NW/) {         # it depends---
    if(!(/Scaffold/)) {         # it depends---
      chomp;
      my @temp = split(/\t/);
      $gene2pos->{$temp[0]} = ($temp[2]+$temp[3])/2;
      if($count!=0) {
        if($temp[1] ne $chr) {
          $vchr->[$index] = $gene2pos->{$gene};
          $index2vchr->{$index-1} = $chr;
          $index += 1;
        }
      }
      else {
        $vchr->[$index] = 0;
        $index += 1;
      }
      $gene2pos->{$temp[0]} += $vchr->[$index-1];

      $count += 1;
      $gene = $temp[0];
      $chr = $temp[1];
    }
  }
  $vchr->[$index] = $gene2pos->{$gene};
  $index2vchr->{$index-1} = $chr;
  $$maxy = $gene2pos->{$gene};

  close(FR1);
  system("rm $out1");
}

sub drawDataF {
  my ($xw,$yh,$x0,$y0,$x1,$y1,$x2,$y2,$im,$maxx,$maxy,$bgene2agene,$gene2pos,$gene2sites,$color) = @_;

  my ($x,$y);
  foreach my $i (keys %$bgene2agene) {
    my ($xx,$yy) = ($gene2sites->{$bgene2agene->{$i}},$gene2sites->{$i});
    my ($xt1,$xt2) = split(/\t/,$xx);
    my ($yt1,$yt2) = split(/\t/,$yy);
    if(abs($yt1-$yt2)>0) {
      my ($xtm,$ytm) = (abs($xt1-$xt2)/2,abs($yt1-$yt2)/2);
      my ($xt0,$yt0) = ($gene2pos->{$bgene2agene->{$i}},$gene2pos->{$i});
      if($xt1<$xt2 && $yt1>$yt2 || $xt1>$xt2 && $yt1<$yt2) {
        ($x1,$x2) = ($xt0-$xtm,$xt0+$xtm);
        ($y1,$y2) = ($yt0+$ytm,$yt0-$ytm);
      }
      else {
        ($x1,$x2) = ($xt0+$xtm,$xt0-$xtm);
        ($y1,$y2) = ($yt0+$ytm,$yt0-$ytm);
      }
      ($x,$y) = ($x1,$y1);
      &transRatio(\$x,\$y,$xw,$yh,$maxx,$maxy);
      ($x1,$y1) = ($x,$y);
      ($x,$y) = ($x2,$y2);
      &transRatio(\$x,\$y,$xw,$yh,$maxx,$maxy);
      ($x2,$y2) = ($x,$y);
      &transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
      #$$im->circle(cx=>$x1,cy=>$y1,r=>0.5,stroke=>$color,'stroke-width'=>0.4,fill=>$color);
      $$im->line(x1=>$x1,y1=>$y1,x2=>$x2,y2=>$y2,stroke=>$color,'stroke-width'=>3,fill=>'red');
      #$$im->circle(cx=>$x2,cy=>$y2,r=>0.5,stroke=>$color,'stroke-width'=>0.4,fill=>$color);
    }
  }
}

sub drawData {
  my ($xw,$yh,$x0,$y0,$x1,$y1,$x2,$y2,$im,$maxx,$maxy,$bgene2agene,$gene2pos,$vchr,$hchr,$color) = @_;

  my ($x,$y);
  foreach my $i (keys %$bgene2agene) {
    #print "11$i\n";
    ($x,$y) = ($gene2pos->{$bgene2agene->{$i}},$gene2pos->{$i});
    &transRatio(\$x,\$y,$xw,$yh,$maxx,$maxy);
    ($x1,$y1) = ($x,$y);
    ($x2,$y2) = ($x,$y);
    &transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
    $$im->circle(cx=>$x1,cy=>$y1,r=>2,stroke=>$color,'stroke-width'=>0.4,fill=>$color);
  }

}

sub extractFrags {
  my ($coord,$gposb,$gposa,$hits) = @_;

  open(my $FR,$coord);
  open(my $FW1,">$hits");
  open(my $FW2,">$gposb");
  open(my $FW3,">$gposa");
  my $index = 0;
  my ($da,$db);
  while($_=<$FR>) {
    s/^\s+//;
    if(/^\d+/) {
      s/\|//g;
      chomp;
      my @temp = split(/\s+/);
      if($temp[11]==1) {
        $da = "+";
      }
      else {
        $da = "-";
      }
      if($temp[12]==1) {
        $db = "+";
      }
      else {
        $db = "-";
      }
      print $FW3 "FragmentA$index\tA$temp[13]\t$temp[0]\t$temp[1]\t$da\n";
      print $FW2 "FragmentB$index\tB$temp[14]\t$temp[2]\t$temp[3]\t$db\n";
      print $FW1 "FragmentB$index\tB$temp[14]\t$temp[2]\t$temp[3]\tFragmentA$index\tA$temp[13]\t$temp[0]\t$temp[1]\t$temp[6]\t*\t*\t*\t*\n";
      $index += 1;
    }
  }
  close($FR);
  close($FW1);
  close($FW2);
  close($FW3);
}

sub processSyn {
  my ($in,$out) = @_;

  open(my $FW,">$out");
  open(my $FR,$in);
  while($_=<$FR>) {
    chomp;
    my @temp = split(/\t/);
    $temp[1] =~s/SL4\.0//;
    $temp[1]=~/^([^\d]+)(\d+)$/;
    my $chr = int($2);
    print $FW "a:$temp[0]\ta:$1";
    if($chr<10) {
      print $FW "00$chr";
    }
    elsif($chr<100) {
      print $FW "0$chr";
    }
    else {
      print $FW "$chr";
    }
    print $FW "\t$temp[2]\t$temp[3]\t";

    $temp[5]=~/^([^\d]+)(\d+)$/;
    $chr = int($2);
    print $FW "b:$temp[4]\tb:$1";
    if($chr<10) {
      print $FW "00$chr";
    }
    elsif($chr<100) {
      print $FW "0$chr";
    }
    else {
      print $FW "$chr";
    }
    print $FW "\t$temp[6]\t$temp[7]\n";
  }
  close($FR);
  close($FW);
}


sub pcolor {
#syntenic points------------------------------
  my %block2color = (
  'A'=>'#F7D471', 'B'=>'#F7D471', 'C'=>'#F7D471',
  'D'=>'#C3253C', 'E'=>'#C3253C',
  'F'=>'#23264F', 'G'=>'#23264F', 'H'=>'#23264F',
  'I'=>'#B55A98', 'J'=>'#B55A98',
  'K'=>'#DC932C', 'L'=>'#DC932C', 'M'=>'#DC932C', 'N'=>'#DC932C',
  'O'=>'#1D923E', 'P'=>'#1D923E', 'Q'=>'#1D923E', 'R'=>'#1D923E',
  'S'=>'#48AACA', 'T'=>'#48AACA', 'U'=>'#48AACA',
  'V'=>'#DA9099', 'W'=>'#DA9099', 'X'=>'#DA9099',
  'Y'=>'white','NULL'=>'black'
  );
  return \%block2color;
}

sub setCanvas {
  my ($widtn,$height,$xw,$yh,$x0,$y0,$im) = @_;
  $$im = SVG->new(width=>$widtn,height=>$height);
  ($$xw,$$yh,$$x0,$$y0) = ($widtn-300,$height-300,200,$height-200);
}

sub drawAxis {
  my ($xw,$yh,$x0,$y0,$x1,$y1,$x2,$y2,$im) = @_;
  my ($x,$y);
  ($x1,$y1,$x2,$y2) = (0,0,$xw,0);
  &transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
  $$im->line(x1=>$x1,y1=>$y1,x2=>$x2,y2=>$y2,stroke=>'black','stroke-width'=>0.5,fill=>'black');
  ($x1,$y1,$x2,$y2) = (0,0,0,$yh);
  &transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
  $$im->line(x1=>$x1,y1=>$y1,x2=>$x2,y2=>$y2,stroke=>'black','stroke-width'=>0.5,fill=>'black');
}

sub drawXYpos {
  my ($xw,$yh,$x0,$y0,$x1,$y1,$x2,$y2,$im,$maxx,$maxy) = @_;
  #axis x----------------
  my ($x,$y);
  ($x1,$y1,$x2,$y2) = (0,0,$xw,0);
  #&transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
  #$im->line(x1=>$x1,y1=>$y1,x2=>$x2,y2=>$y2,stroke=>'black','stroke-width'=>1,fill=>'black');
  my $flag = 1;
  for(my $i=30;$i*1000000<$maxx;$i+=30) {
    ($x,$y) = ($i*1000000,0);
    &transRatio(\$x,\$y,$xw,$yh,$maxx,$maxy);
    ($x1,$y1,$x2,$y2) = ($x,0,$x,-10);
    &transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
    $$im->line(x1=>$x1,y1=>$y1,x2=>$x2,y2=>$y2,stroke=>'black','stroke-width'=>1,fill=>'black');
    if($flag==1) {
      ($x1,$y1) = ($x-20,-10-15);
      &transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
      $$im->text(
        x=>$x1, y=>$y1,
        style => {
        'font-family' => 'Arial',
        'font-size' => '14pt',
        #'stroke' => 'black',
        'font-weight' => 'bold'
        }
      )->cdata($i." Mb");
    }
    if($flag==1) {
      $flag = 0;
    }
    else {
      $flag = 1;
    }
  }

  #axis y----------------
  ($x1,$y1,$x2,$y2) = (0,0,0,$yh);
  #&transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
  #$im->line(x1=>$x1,y1=>$y1,x2=>$x2,y2=>$y2,stroke=>'black','stroke-width'=>1,fill=>'black');
  $flag = 1;
  for(my $i=50;$i*1000000<$maxy;$i+=50) {
    ($x,$y) = (0,$i*1000000);
    &transRatio(\$x,\$y,$xw,$yh,$maxx,$maxy);
    ($x1,$y1,$x2,$y2) = (-10,$y,0,$y);
    &transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
    $$im->line(x1=>$x1,y1=>$y1,x2=>$x2,y2=>$y2,stroke=>'black','stroke-width'=>1,fill=>'black');
    if($flag==1) {
      ($x1,$y1) = (-10-8,$y-15-8);
      &transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
      $$im->text(
        x=>$x1, y=>$y1,
        style => {
        'font-family'  => 'Arial',
        #'stroke' => 'black',
        'font-size' => '14pt',
        'font-weight' => 'bold'
        },
        transform => "rotate(-90,$x1,$y1)"
      )->cdata($i." Mb");
    }
    if($flag==1) {
      $flag = 0;
    }
    else {
      $flag = 1;
    }
  }

  ($x1,$y1,$x2,$y2) = (0,$yh,$xw,$yh);
  &transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
  $$im->line(x1=>$x1,y1=>$y1,x2=>$x2,y2=>$y2,stroke=>'black','stroke-width'=>1,fill=>'black');
  ($x1,$y1,$x2,$y2) = ($xw,0,$xw,$yh);
  &transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
  $$im->line(x1=>$x1,y1=>$y1,x2=>$x2,y2=>$y2,stroke=>'black','stroke-width'=>1,fill=>'black');
  #$im->rect(x1=>$x1x,y1=>$y1x,height=>1000,width=>2000,fill=>'threedlightshadow');
}

sub csep {
  my ($xw,$yh,$x0,$y0,$x1,$y1,$x2,$y2,$im,$maxx,$maxy,$vchr,$hchr,$index2vchr,$index2hchr) = @_;
  my ($x,$y);

  #chromosome sepatater of At-------------------
  my $i;
  for($i=1;$i<@$vchr;$i+=1) {
    ($x,$y) = (0,$vchr->[$i]);
    &transRatio(\$x,\$y,$xw,$yh,$maxx,$maxy);
    ($x1,$y1,$x2,$y2) = (0,$y,$xw,$y);
    &transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
    $$im->line(x1=>$x1,y1=>$y1,x2=>$x2,y2=>$y2,stroke=>'black','stroke-width'=>0.2,fill=>'black');

    ($x,$y) = (0,($vchr->[$i]+$vchr->[$i-1])/2);
    &transRatio(\$x,\$y,$xw,$yh,$maxx,$maxy);
    ($x1,$y1) = (-10-30-80,$y-15-35);
    &transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
    $$im->text(
      x=>$x1, y=>$y1,
      style => {
      'font-family'  => 'Arial',
      #'stroke' => 'black',
      'font-size' => '14pt',
      'font-weight' => 'bold'
      },
      #transform => "rotate(-90,$x1,$y1)"
      transform => "rotate(-45,$x1,$y1)"
    )->cdata($index2vchr->{$i-1});
  }
  
  #chromosome sepatater of Br-------------------
  for($i=1;$i<@$hchr;$i+=1) {
    ($x,$y) = ($hchr->[$i],0);
    &transRatio(\$x,\$y,$xw,$yh,$maxx,$maxy);
    ($x1,$y1,$x2,$y2) = ($x,0,$x,$yh);
    &transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
    $$im->line(x1=>$x1,y1=>$y1,x2=>$x2,y2=>$y2,stroke=>'black','stroke-width'=>1,fill=>'black');

    ($x,$y) = (($hchr->[$i]+$hchr->[$i-1])/2,0);
    &transRatio(\$x,\$y,$xw,$yh,$maxx,$maxy);
    ($x1,$y1) = ($x-70,-10-45-5-50-50);
    &transPos(\$x1,\$y1,\$x2,\$y2,$x0,$y0);
    $$im->text(
      x=>$x1, y=>$y1,
      style => {
      'font-family'  => 'Arial',
      #'stroke' => 'black',
      'font-size' => '14pt',
      'font-weight' => 'bold'
      },
      transform => "rotate(-45,$x1,$y1)"
    )->cdata($index2hchr->{$i-1});
  }
}

sub transPos {
  my ($x1,$y1,$x2,$y2,$x0,$y0) = @_;
  $$x1 += $x0;
  $$y1 = $y0 - $$y1;
  $$x2 += $x0;
  $$y2 = $y0 - $$y2;
}

sub transRatio {
  my ($x,$y,$xw,$yh,$maxx,$maxy) = @_;
  #print "$$x\t$maxx\n";
  $$x /= $maxx;
  $$x *= $xw;
  $$y /= $maxy;
  $$y *= $yh;
}

