|
0
|
1 /*****************************************************************************
|
|
|
2 flankBed.cpp
|
|
|
3
|
|
|
4 (c) 2009 - Aaron Quinlan
|
|
|
5 Hall Laboratory
|
|
|
6 Department of Biochemistry and Molecular Genetics
|
|
|
7 University of Virginia
|
|
|
8 aaronquinlan@gmail.com
|
|
|
9
|
|
|
10 Licensed under the GNU General Public License 2.0 license.
|
|
|
11 ******************************************************************************/
|
|
|
12 #include "lineFileUtilities.h"
|
|
|
13 #include "flankBed.h"
|
|
|
14
|
|
|
15
|
|
|
16 BedFlank::BedFlank(string &bedFile, string &genomeFile, bool forceStrand, float leftFlank, float rightFlank, bool fractional) {
|
|
|
17
|
|
|
18 _bedFile = bedFile;
|
|
|
19 _genomeFile = genomeFile;
|
|
|
20 _forceStrand = forceStrand;
|
|
|
21 _leftFlank = leftFlank;
|
|
|
22 _rightFlank = rightFlank;
|
|
|
23 _fractional = fractional;
|
|
|
24
|
|
|
25 _bed = new BedFile(bedFile);
|
|
|
26 _genome = new GenomeFile(genomeFile);
|
|
|
27
|
|
|
28 // get going, slop it up.
|
|
|
29 FlankBed();
|
|
|
30 }
|
|
|
31
|
|
|
32
|
|
|
33 BedFlank::~BedFlank(void) {
|
|
|
34
|
|
|
35 }
|
|
|
36
|
|
|
37
|
|
|
38 void BedFlank::FlankBed() {
|
|
|
39
|
|
|
40 int lineNum = 0;
|
|
|
41 BED bedEntry, nullBed; // used to store the current BED line from the BED file.
|
|
|
42 BedLineStatus bedStatus;
|
|
|
43
|
|
|
44 _bed->Open();
|
|
|
45 bedStatus = _bed->GetNextBed(bedEntry, lineNum);
|
|
|
46 while (bedStatus != BED_INVALID) {
|
|
|
47 if (bedStatus == BED_VALID) {
|
|
|
48
|
|
|
49 int leftFlank = _leftFlank;
|
|
|
50 int rightFlank = _rightFlank;
|
|
|
51 if (_fractional == true) {
|
|
|
52 leftFlank = (int) (_leftFlank * bedEntry.size());
|
|
|
53 rightFlank = (int) (_rightFlank * bedEntry.size());
|
|
|
54 }
|
|
|
55
|
|
|
56 if ((_forceStrand == false) || (bedEntry.strand == "+"))
|
|
|
57 {
|
|
|
58 AddFlank(bedEntry, leftFlank, rightFlank);
|
|
|
59 }
|
|
|
60 else if ((_forceStrand == true) && (bedEntry.strand == "-" ))
|
|
|
61 {
|
|
|
62 AddStrandedFlank(bedEntry, leftFlank, rightFlank);
|
|
|
63 }
|
|
|
64 bedEntry = nullBed;
|
|
|
65 }
|
|
|
66 bedStatus = _bed->GetNextBed(bedEntry, lineNum);
|
|
|
67 }
|
|
|
68 _bed->Close();
|
|
|
69 }
|
|
|
70
|
|
|
71
|
|
|
72 void BedFlank::AddFlank(BED &bed, int leftFlank, int rightFlank) {
|
|
|
73
|
|
|
74 int chromSize = _genome->getChromSize(bed.chrom);
|
|
|
75 if (chromSize == -1) {
|
|
|
76 cerr << "ERROR: chrom \"" << bed.chrom << "\" not found in genome file. Exiting." << endl;
|
|
|
77 exit(1);
|
|
|
78 }
|
|
|
79
|
|
|
80 // init. our left and right flanks to the original BED entry.
|
|
|
81 // we'll create the flanks from these coordinates.
|
|
|
82 BED left = bed;
|
|
|
83 BED right = bed;
|
|
|
84
|
|
|
85 // make the left flank (if necessary)
|
|
|
86 if (leftFlank > 0) {
|
|
|
87 if ( (static_cast<int>(left.start) - leftFlank) > 0)
|
|
|
88 {
|
|
|
89 left.end = left.start;
|
|
|
90 left.start -= leftFlank;
|
|
|
91 }
|
|
|
92 else
|
|
|
93 {
|
|
|
94 left.end = left.start;
|
|
|
95 left.start = 0;
|
|
|
96 }
|
|
|
97 // report the left flank
|
|
|
98 _bed->reportBedNewLine(left);
|
|
|
99 }
|
|
|
100
|
|
|
101 // make the left flank (if necessary)
|
|
|
102 if (rightFlank > 0) {
|
|
|
103 if ( (static_cast<int>(right.end) + (rightFlank+1)) <= static_cast<int>(chromSize))
|
|
|
104 {
|
|
|
105 right.start = right.end;
|
|
|
106 right.end += (rightFlank);
|
|
|
107 }
|
|
|
108 else {
|
|
|
109 right.start = right.end;
|
|
|
110 right.end += chromSize;
|
|
|
111 }
|
|
|
112 // report the right flank
|
|
|
113 _bed->reportBedNewLine(right);
|
|
|
114 }
|
|
|
115 }
|
|
|
116
|
|
|
117
|
|
|
118 void BedFlank::AddStrandedFlank(BED &bed, int leftFlank, int rightFlank) {
|
|
|
119
|
|
|
120 int chromSize = _genome->getChromSize(bed.chrom);
|
|
|
121 if (chromSize == -1) {
|
|
|
122 cerr << "ERROR: chrom \"" << bed.chrom << "\" not found in genome file. Exiting." << endl;
|
|
|
123 exit(1);
|
|
|
124 }
|
|
|
125
|
|
|
126 // init. our left and right flanks to the original BED entry.
|
|
|
127 // we'll create the flanks from these coordinates.
|
|
|
128 BED left = bed;
|
|
|
129 BED right = bed;
|
|
|
130
|
|
|
131 // make the left flank (if necessary)
|
|
|
132 if (rightFlank > 0) {
|
|
|
133 if ( (static_cast<int>(left.start) - rightFlank) > 0)
|
|
|
134 {
|
|
|
135 left.end = left.start;
|
|
|
136 left.start -= rightFlank;
|
|
|
137 }
|
|
|
138 else
|
|
|
139 {
|
|
|
140 left.end = left.start;
|
|
|
141 left.start = 0;
|
|
|
142 }
|
|
|
143 // report the left flank
|
|
|
144 _bed->reportBedNewLine(left);
|
|
|
145 }
|
|
|
146
|
|
|
147 // make the left flank (if necessary)
|
|
|
148 if (leftFlank > 0) {
|
|
|
149 if ( (static_cast<int>(right.end) + leftFlank) <= static_cast<int>(chromSize))
|
|
|
150 {
|
|
|
151 right.start = right.end;
|
|
|
152 right.end += leftFlank;
|
|
|
153 }
|
|
|
154 else {
|
|
|
155 right.start = right.end;
|
|
|
156 right.end = chromSize;
|
|
|
157 }
|
|
|
158 // report the right flank
|
|
|
159 _bed->reportBedNewLine(right);
|
|
|
160 }
|
|
|
161 }
|
|
|
162
|
|
|
163
|
