diff scripts/make.bp.bed.pl @ 0:28d1a6f8143f draft

planemo upload for repository https://github.com/portiahollyoak/Tools commit 132bb96bba8e7aed66a102ed93b7744f36d10d37-dirty

--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/scripts/make.bp.bed.pl	Mon Apr 25 13:08:56 2016 -0400
@@ -0,0 +1,110 @@
+#! /usr/bin/perl
+
+use strict;
+
+my @sample=();
+open (in, "<$ARGV[0]") or die "Can't open $ARGV[0] since $!\n";
+my $line=<in>;
+close in;
+
+my %chrs=();
+my @a=split(/\t/, $line);
+for my $i (0..$#a) {
+    if ($a[$i] =~ /_class$/) {
+	my $name=$a[$i];
+	$name =~ s/_class//;
+	my $j=$i+1;
+	my $k=$i+2;
+	my $l=$i+3;
+	system("cut -f7,4,6,$j,$k,$l $ARGV[0] > temp");
+	open (input, "<temp") or die "Can't open temp since $!\n";
+	open (output, ">>$name.insertion.bp.bed") or die "Can't open $name.insertion.bp.bed since $!\n";
+	my $header=<input>;
+	while (my $line=<input>) {
+	    chomp($line);
+	    my @b=split(/\t/, $line);
+	    if (($b[4] ne "0")||($b[5] ne "0")) {
+		my @c=split(/\:/, $b[2]);
+		my @d=split(/\./, $c[1]);
+		if ($d[0] > $d[1]) {
+		    my $temp=$d[0];
+		    $d[0]=$d[1];
+		    $d[1]=$temp;
+		}
+		my $lower=$d[0];
+		my $upper=$d[1];
+	        if (($lower >= 0) && ($upper >= 0)) {
+		   print output "$c[0]\t$lower\t$upper\t$b[0]\t$b[1]\t$b[3]\t$b[4]\t$b[5]\n";
+	        }
+		$chrs{$c[0]}=1;
+	    }
+	}
+	close input;
+	close output;
+	system("rm temp");
+
+	if ($ARGV[1] ne "") {
+	    open (input, "<$name.insertion.bp.bed") or die "Can't open $name.insertion.bp.bed since $!\n";
+	    open (output, ">tmp") or die "Can't tmp since $!\n";
+	    while (my $line=<input>) {
+		chomp($line);
+		my @a=split(/\t/, $line);
+		if (($a[0] =~ /^\d{1,2}$/) || ($a[0] eq "X") || ($a[0] eq "Y")) {$a[0]="chr$a[0]";}
+		my $strand="+";
+		if ($a[4] eq "antisense") {$strand="-";}
+		print output "$a[0]\t$a[1]\t$a[2]\t$a[3]\t\.\t$strand\t$a[5]\t$a[6]\t$a[7]\n";
+	    }
+	    close input;
+	    close output;
+
+	    system("bedtools intersect -a tmp -b $ARGV[1] -f 0.1 -wo -s > tmp1");
+	    if ($ARGV[2] eq "") {
+		system("awk -F \"\\t\" '{OFS=\"\\t\"; if ((\$4==\$13)&&(\$6==\$15)) print \$1,\$2,\$3,\$4,\$5,\$6}' tmp1 > tmp2");
+	    }
+	    else {
+		my %family=();
+		open (input, "<$ARGV[2]") or die "Can't open $ARGV[2] since $!\n";
+		while (my $line=<input>) {
+		    chomp($line);
+		    my @a=split(/\t/, $line);
+		    $family{$a[0]}=$a[1];
+		}
+		close input;
+
+		open (input, "<tmp1") or die "Can't open tmp1 since $!\n";
+		open (output, ">>tmp2") or die "Can't open tmp2 since $!\n";
+		while (my $line=<input>) {
+		    chomp($line);
+		    my @a=split(/\t/, $line);
+		    if (($family{$a[3]} eq $family{$a[12]}) && ($a[5] eq $a[14])) {
+			print output "$a[0]\t$a[1]\t$a[2]\t$a[3]\t$a[4]\t$a[5]\n";
+		    }
+		}
+		close input;
+		close output;
+	    }
+	    
+	    if (-s "tmp2") {
+		system("bedtools subtract -a tmp -b tmp2 -f 1.0 > tmp3");
+		open (input, "<tmp3") or die "Can't open tmp3 since $!\n";
+		open (output, ">$name.insertion.bp.bed") or die "Can't open $name.insertion.bp.bed since $!\n";
+		while (my $line=<input>) {
+		    chomp($line);
+		    my @a=split(/\t/, $line);
+		    my $direction="sense";
+		    if ($a[5] eq "-") {$direction="antisense";}
+		    my $chr_num=$a[0];
+		    $chr_num =~ s/chr//;
+		    if (($chrs{$a[0]} == 1) && (! defined $chrs{$chr_num})) {$chr_num=$a[0];}
+		    print output "$chr_num\t$a[1]\t$a[2]\t$a[3]\t$direction\t$a[6]\t$a[7]\t$a[8]\n";
+		}
+		close input;
+		close output;
+	    }
+	}
+
+	system("rm tmp*");
+
+    }
+}
+