File: prop_pair_sam_refiner.pl

package info (click to toggle)
trinityrnaseq 2.11.0%2Bdfsg-6
  • links: PTS, VCS
  • area: main
  • in suites: bullseye
  • size: 417,528 kB
  • sloc: perl: 48,420; cpp: 17,749; java: 12,695; python: 3,124; sh: 1,030; ansic: 983; makefile: 688; xml: 62
file content (102 lines) | stat: -rwxr-xr-x 1,972 bytes parent folder | download | duplicates (2)
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
#!/usr/bin/env perl

use strict;
use warnings;

use lib ("/usr/lib/trinityrnaseq/PerlLib");
use SAM_reader;
use SAM_entry;

my $usage = "\n\nusage: name_sorted_paired_reads.sam\n\n";

my $sam_file = $ARGV[0] or die $usage;

main: {


    my $prev_read_name = "";
    my $prev_scaff_name = "";

    my @reads;

    my $sam_reader = new SAM_reader($sam_file);
    while ($sam_reader->has_next()) {
        

        my $read = $sam_reader->get_next();
        
        my $scaff_name = $read->get_scaffold_name();
        my $core_read_name = $read->get_core_read_name();
        
        if ($scaff_name ne $prev_scaff_name ||  $core_read_name ne $prev_read_name) {

            if (@reads) {
                &process_pairs(@reads);
                @reads = ();
            }
        }
        
        push (@reads, $read);
        

        $prev_read_name = $core_read_name;
        $prev_scaff_name = $scaff_name;
        
        
    }
    
    &process_pairs(@reads);


    exit(0);
}

####
sub process_pairs {
    my (@reads) = @_;


    my @left_reads;
    my @right_reads;

    foreach my $read (@reads) {
        if ($read->is_first_in_pair()) {
            push (@left_reads, $read);
        }
        elsif ($read->is_second_in_pair()) {
            push (@right_reads, $read);
        }
    }

    unless (@left_reads && @right_reads) {
        die "Error, dont have pairs!";
    }
    
    my $left_read = shift @left_reads;
    my $aligned_pos = $left_read->get_aligned_position();

    my $right_read = undef;
    while ((! defined($right_read)) && @right_reads) {
        my $read = shift @right_reads;
        if ($read->get_mate_scaffold_position() == $aligned_pos) {
            $right_read = $read;
            last;
        }
    }

    unless ($left_read && $right_read) {
        die "Error, couldn't match pairs.";
    }

    print $left_read->toString() . "\n";
    print $right_read->toString() . "\n";
    

    return;
}