-
Notifications
You must be signed in to change notification settings - Fork 0
/
takeinteractionblock.pl
executable file
·104 lines (79 loc) · 2.12 KB
/
takeinteractionblock.pl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
#!/usr/bin/env perl
use strict;
use warnings;
# Warning: this script, as it is, may use huges amount of RAM
# @l arrays may be tied to temporary files instead of packed strings
BEGIN {
use FindBin '$Bin';
use Tie::Array::Packed;
require "$Bin/share/Metaheader.pm";
}
if ($#ARGV != 1) {
print STDERR "usage: ./takeinteractionblock.pl block1 block2 "
, "< input > output\n";
exit -1;
}
my $blockstring1 = shift @ARGV;
my $blockstring2 = shift @ARGV;
# Read the header
my $header = <>;
chomp($header);
my $metah = Metaheader->new($header);
my @rowarray = @{ $metah->{rowinfo} };
my @output1 = $metah->selectvector($blockstring1);
my @output2 = $metah->selectvector($blockstring2);
die "No match for block 1" if ($#output1 == -1);
die "No match for block 2" if ($#output2 == -1);
# check who is first
my $unofirst = $output1[0] < $output2[0] ? 1 : 0;
# check for overlap
die "Blocks are overlapping" if
( $unofirst ?
$output2[0] < $output1[$#output1] :
$output1[0] < $output2[$#output2] );
if ($unofirst) {
# print header
print join("\t", "\"die\"",
@{ $metah->{strings} }[@output2]), "\n";
my $j = 0;
while (<>) {
last if ($#output1 == -1);
if ($j == $output1[0]) {
chomp;
my @r=split("\t");
my $head=shift @r;
my @rsel = map { $_ - $j } @output2;
print join("\t", $head, @r[@rsel]),"\n";
shift @output1;
}
$j++;
}
} else {
# not $unofirst
# print header
print join("\t", "\"die\"",
@{ $metah->{strings} }[@output2]), "\n";
# Warning: this script, as it is, may use huges amount of RAM
# @l arrays may be tied to temporary files
my @l = map {Tie::Array::Packed::IntegerNative->make()} @output2;
my $j = 0;
while (<>) {
last if ($#output2 == -1);
if ($j == $output2[0]) {
chomp;
my @r=split("\t");
my $head=shift @r;
my @rsel = map { $_ - $j } @output1;
for my $i (0..$#output1) {
push @{$l[$i]}, $r[$output1[$i]-$j];
}
shift @output2;
}
$j++;
}
for my $i (0..$#output1) {
print join("\t", $metah->{strings}->[$output1[$i]],
@{ $l[$i] }),"\n";
}
}
0;