libStatGen Software  1
TestCigarHelper.cpp
1 /*
2  * Copyright (C) 2011 Regents of the University of Michigan
3  *
4  * This program is free software: you can redistribute it and/or modify
5  * it under the terms of the GNU General Public License as published by
6  * the Free Software Foundation, either version 3 of the License, or
7  * (at your option) any later version.
8  *
9  * This program is distributed in the hope that it will be useful,
10  * but WITHOUT ANY WARRANTY; without even the implied warranty of
11  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12  * GNU General Public License for more details.
13  *
14  * You should have received a copy of the GNU General Public License
15  * along with this program. If not, see <http://www.gnu.org/licenses/>.
16  */
17 
18 #include "TestCigarHelper.h"
19 #include "TestValidate.h"
20 #include "CigarHelper.h"
21 #include <assert.h>
22 
23 void testCigarHelper()
24 {
25  // Call generic test.
26  CigarHelperTest::testCigarHelper();
27 }
28 
29 
30 void CigarHelperTest::testCigarHelper()
31 {
32  testSoftClipBeginByRefPos();
33  testSoftClipEndByRefPos();
34 }
35 
36 
37 void CigarHelperTest::testSoftClipBeginByRefPos()
38 {
39  SamRecord record;
40  CigarRoller newCigar;
41  std::string newCigarString;
42  int32_t newPos = 0;
43 
44  // Setup the current Cigar.
45  // Cigar: HHHSSSMMMDDDMMMIIIMMMPPPMMMDDDMMMSSSHHH
46  // ReadPos: 000000 000011111 111 112222
47  // ReadPos: 012345 678901234 567 890123
48  // RefPos: 111111111 122 222222223
49  // RefPos: 012345678 901 234567890
50  const char* origCigar = "3H3S3M3D3M3I3M3P3M3D3M3S3H";
51  record.setCigar(origCigar);
52  record.set0BasedPosition(10);
53  record.setSequence("gggAAATTTCCCTTTGGGAAAggg");
54 
55  ////////////////////////////////////////////////////////
56  // Clip outside of the range (after). Everything should be clipped.
57  assert(CigarHelper::softClipBeginByRefPos(record, 10000, newCigar, newPos) == 23);
58  newCigar.getCigarString(newCigarString);
59  assert(strcmp(newCigarString.c_str(), "3H24S3H") == 0);
60 
61  ////////////////////////////////////////////////////////
62  // Clip outside of the range (before). Nothing should change.
63  assert(CigarHelper::softClipBeginByRefPos(record, 1, newCigar, newPos) ==
64  CigarHelper::NO_CLIP);
65  newCigar.getCigarString(newCigarString);
66  assert(strcmp(newCigarString.c_str(), origCigar) == 0);
67 
68 
69  ////////////////////////////////////////////////////////
70  ////////////////////////////////////////////////////////
71  // Test clipping at every position of the read.
72 
73  ////////////////////////////////////////////////////////
74  // Clip at the first position.
75  assert(CigarHelper::softClipBeginByRefPos(record, 10, newCigar, newPos) == 3);
76  assert(newPos == 11);
77  newCigar.getCigarString(newCigarString);
78  //std::cout << newCigarString.c_str() << std::endl;
79  assert(strcmp(newCigarString.c_str(), "3H4S2M3D3M3I3M3P3M3D3M3S3H") == 0);
80 
81  ////////////////////////////////////////////////////////
82  // Clip in the middle of the first Match.
83  assert(CigarHelper::softClipBeginByRefPos(record, 11, newCigar, newPos) == 4);
84  assert(newPos == 12);
85  newCigar.getCigarString(newCigarString);
86  //std::cout << newCigarString.c_str() << std::endl;
87  assert(strcmp(newCigarString.c_str(), "3H5S1M3D3M3I3M3P3M3D3M3S3H") == 0);
88 
89  ////////////////////////////////////////////////////////
90  assert(CigarHelper::softClipBeginByRefPos(record, 12, newCigar, newPos) == 5);
91  assert(newPos == 16);
92  newCigar.getCigarString(newCigarString);
93  //std::cout << newCigarString.c_str() << std::endl;
94  assert(strcmp(newCigarString.c_str(), "3H6S3M3I3M3P3M3D3M3S3H") == 0);
95 
96  ////////////////////////////////////////////////////////
97  assert(CigarHelper::softClipBeginByRefPos(record, 13, newCigar, newPos) == 5);
98  assert(newPos == 16);
99  newCigar.getCigarString(newCigarString);
100  //std::cout << newCigarString.c_str() << std::endl;
101  assert(strcmp(newCigarString.c_str(), "3H6S3M3I3M3P3M3D3M3S3H") == 0);
102 
103  ////////////////////////////////////////////////////////
104  assert(CigarHelper::softClipBeginByRefPos(record, 14, newCigar, newPos) == 5);
105  assert(newPos == 16);
106  newCigar.getCigarString(newCigarString);
107  //std::cout << newCigarString.c_str() << std::endl;
108  assert(strcmp(newCigarString.c_str(), "3H6S3M3I3M3P3M3D3M3S3H") == 0);
109 
110  ////////////////////////////////////////////////////////
111  assert(CigarHelper::softClipBeginByRefPos(record, 15, newCigar, newPos) == 5);
112  assert(newPos == 16);
113  newCigar.getCigarString(newCigarString);
114  //std::cout << newCigarString.c_str() << std::endl;
115  assert(strcmp(newCigarString.c_str(), "3H6S3M3I3M3P3M3D3M3S3H") == 0);
116 
117  ////////////////////////////////////////////////////////
118  assert(CigarHelper::softClipBeginByRefPos(record, 16, newCigar, newPos) == 6);
119  assert(newPos == 17);
120  newCigar.getCigarString(newCigarString);
121  //std::cout << newCigarString.c_str() << std::endl;
122  assert(strcmp(newCigarString.c_str(), "3H7S2M3I3M3P3M3D3M3S3H") == 0);
123 
124  ////////////////////////////////////////////////////////
125  assert(CigarHelper::softClipBeginByRefPos(record, 17, newCigar, newPos) == 7);
126  assert(newPos == 18);
127  newCigar.getCigarString(newCigarString);
128  //std::cout << newCigarString.c_str() << std::endl;
129  assert(strcmp(newCigarString.c_str(), "3H8S1M3I3M3P3M3D3M3S3H") == 0);
130 
131  ////////////////////////////////////////////////////////
132  assert(CigarHelper::softClipBeginByRefPos(record, 18, newCigar, newPos) == 11);
133  assert(newPos == 19);
134  newCigar.getCigarString(newCigarString);
135  //std::cout << newCigarString.c_str() << std::endl;
136  assert(strcmp(newCigarString.c_str(), "3H12S3M3P3M3D3M3S3H") == 0);
137 
138  ////////////////////////////////////////////////////////
139  assert(CigarHelper::softClipBeginByRefPos(record, 19, newCigar, newPos) == 12);
140  assert(newPos == 20);
141  newCigar.getCigarString(newCigarString);
142  //std::cout << newCigarString.c_str() << std::endl;
143  assert(strcmp(newCigarString.c_str(), "3H13S2M3P3M3D3M3S3H") == 0);
144 
145  ////////////////////////////////////////////////////////
146  assert(CigarHelper::softClipBeginByRefPos(record, 20, newCigar, newPos) == 13);
147  assert(newPos == 21);
148  newCigar.getCigarString(newCigarString);
149  //std::cout << newCigarString.c_str() << std::endl;
150  assert(strcmp(newCigarString.c_str(), "3H14S1M3P3M3D3M3S3H") == 0);
151 
152  ////////////////////////////////////////////////////////
153  assert(CigarHelper::softClipBeginByRefPos(record, 21, newCigar, newPos) == 14);
154  assert(newPos == 22);
155  newCigar.getCigarString(newCigarString);
156  //std::cout << newCigarString.c_str() << std::endl;
157  assert(strcmp(newCigarString.c_str(), "3H15S3M3D3M3S3H") == 0);
158 
159  ////////////////////////////////////////////////////////
160  assert(CigarHelper::softClipBeginByRefPos(record, 22, newCigar, newPos) == 15);
161  assert(newPos == 23);
162  newCigar.getCigarString(newCigarString);
163  //std::cout << newCigarString.c_str() << std::endl;
164  assert(strcmp(newCigarString.c_str(), "3H16S2M3D3M3S3H") == 0);
165 
166  ////////////////////////////////////////////////////////
167  assert(CigarHelper::softClipBeginByRefPos(record, 23, newCigar, newPos) == 16);
168  assert(newPos == 24);
169  newCigar.getCigarString(newCigarString);
170  //std::cout << newCigarString.c_str() << std::endl;
171  assert(strcmp(newCigarString.c_str(), "3H17S1M3D3M3S3H") == 0);
172 
173  ////////////////////////////////////////////////////////
174  assert(CigarHelper::softClipBeginByRefPos(record, 24, newCigar, newPos) == 17);
175  assert(newPos == 28);
176  newCigar.getCigarString(newCigarString);
177  //std::cout << newCigarString.c_str() << std::endl;
178  assert(strcmp(newCigarString.c_str(), "3H18S3M3S3H") == 0);
179 
180  ////////////////////////////////////////////////////////
181  assert(CigarHelper::softClipBeginByRefPos(record, 25, newCigar, newPos) == 17);
182  assert(newPos == 28);
183  newCigar.getCigarString(newCigarString);
184  //std::cout << newCigarString.c_str() << std::endl;
185  assert(strcmp(newCigarString.c_str(), "3H18S3M3S3H") == 0);
186 
187  ////////////////////////////////////////////////////////
188  assert(CigarHelper::softClipBeginByRefPos(record, 26, newCigar, newPos) == 17);
189  assert(newPos == 28);
190  newCigar.getCigarString(newCigarString);
191  //std::cout << newCigarString.c_str() << std::endl;
192  assert(strcmp(newCigarString.c_str(), "3H18S3M3S3H") == 0);
193 
194  ////////////////////////////////////////////////////////
195  assert(CigarHelper::softClipBeginByRefPos(record, 27, newCigar, newPos) == 17);
196  assert(newPos == 28);
197  newCigar.getCigarString(newCigarString);
198  //std::cout << newCigarString.c_str() << std::endl;
199  assert(strcmp(newCigarString.c_str(), "3H18S3M3S3H") == 0);
200 
201  ////////////////////////////////////////////////////////
202  assert(CigarHelper::softClipBeginByRefPos(record, 28, newCigar, newPos) == 18);
203  assert(newPos == 29);
204  newCigar.getCigarString(newCigarString);
205  //std::cout << newCigarString.c_str() << std::endl;
206  assert(strcmp(newCigarString.c_str(), "3H19S2M3S3H") == 0);
207 
208  ////////////////////////////////////////////////////////
209  assert(CigarHelper::softClipBeginByRefPos(record, 29, newCigar, newPos) == 19);
210  assert(newPos == 30);
211  newCigar.getCigarString(newCigarString);
212  //std::cout << newCigarString.c_str() << std::endl;
213  assert(strcmp(newCigarString.c_str(), "3H20S1M3S3H") == 0);
214 
215  ////////////////////////////////////////////////////////
216  assert(CigarHelper::softClipBeginByRefPos(record, 30, newCigar, newPos) == 23);
217  assert(newPos == 10);
218  newCigar.getCigarString(newCigarString);
219  //std::cout << newCigarString.c_str() << std::endl;
220  assert(strcmp(newCigarString.c_str(), "3H24S3H") == 0);
221 
222  ////////////////////////////////////////////////////////
223  assert(CigarHelper::softClipBeginByRefPos(record, 31, newCigar, newPos) == 23);
224  assert(newPos == 10);
225  newCigar.getCigarString(newCigarString);
226  //std::cout << newCigarString.c_str() << std::endl;
227  assert(strcmp(newCigarString.c_str(), "3H24S3H") == 0);
228 
229  ////////////////////////////////////////////////////////
230  ////////////////////////////////////////////////////////
231  // Test clipping at every position when insertions & deletions
232  // are next to each other.
233  origCigar = "3M3D3I3M";
234  record.setCigar(origCigar);
235  record.setSequence("GGGAAAGGG");
236  // Cigar: MMMDDDIIIMMM
237  // ReadPos: 000 000000
238  // ReadPos: 012 345678
239  // RefPos: 111111 111
240  // RefPos: 012345 678
241  record.setCigar(origCigar);
242  assert(CigarHelper::softClipBeginByRefPos(record, 9, newCigar, newPos) ==
243  CigarHelper::NO_CLIP);
244  assert(newPos == 10);
245  newCigar.getCigarString(newCigarString);
246  //std::cout << newCigarString.c_str() << std::endl;
247  assert(strcmp(newCigarString.c_str(), origCigar) == 0);
248 
249  record.setCigar(origCigar);
250  assert(CigarHelper::softClipBeginByRefPos(record, 10, newCigar, newPos) == 0);
251  assert(newPos == 11);
252  newCigar.getCigarString(newCigarString);
253  //std::cout << newCigarString.c_str() << std::endl;
254  assert(strcmp(newCigarString.c_str(), "1S2M3D3I3M") == 0);
255 
256  record.setCigar(origCigar);
257  assert(CigarHelper::softClipBeginByRefPos(record, 11, newCigar, newPos) == 1);
258  assert(newPos == 12);
259  newCigar.getCigarString(newCigarString);
260  //std::cout << newCigarString.c_str() << std::endl;
261  assert(strcmp(newCigarString.c_str(), "2S1M3D3I3M") == 0);
262 
263  record.setCigar(origCigar);
264  assert(CigarHelper::softClipBeginByRefPos(record, 12, newCigar, newPos) == 5);
265  assert(newPos == 16);
266  newCigar.getCigarString(newCigarString);
267  //std::cout << newCigarString.c_str() << std::endl;
268  assert(strcmp(newCigarString.c_str(), "6S3M") == 0);
269 
270  record.setCigar(origCigar);
271  assert(CigarHelper::softClipBeginByRefPos(record, 13, newCigar, newPos) == 5);
272  assert(newPos == 16);
273  newCigar.getCigarString(newCigarString);
274  //std::cout << newCigarString.c_str() << std::endl;
275  assert(strcmp(newCigarString.c_str(), "6S3M") == 0);
276 
277  record.setCigar(origCigar);
278  assert(CigarHelper::softClipBeginByRefPos(record, 14, newCigar, newPos) == 5);
279  assert(newPos == 16);
280  newCigar.getCigarString(newCigarString);
281  //std::cout << newCigarString.c_str() << std::endl;
282  assert(strcmp(newCigarString.c_str(), "6S3M") == 0);
283 
284  record.setCigar(origCigar);
285  assert(CigarHelper::softClipBeginByRefPos(record, 15, newCigar, newPos) == 5);
286  assert(newPos == 16);
287  newCigar.getCigarString(newCigarString);
288  //std::cout << newCigarString.c_str() << std::endl;
289  assert(strcmp(newCigarString.c_str(), "6S3M") == 0);
290 
291  record.setCigar(origCigar);
292  assert(CigarHelper::softClipBeginByRefPos(record, 16, newCigar, newPos) == 6);
293  assert(newPos == 17);
294  newCigar.getCigarString(newCigarString);
295  //std::cout << newCigarString.c_str() << std::endl;
296  assert(strcmp(newCigarString.c_str(), "7S2M") == 0);
297 
298  record.setCigar(origCigar);
299  assert(CigarHelper::softClipBeginByRefPos(record, 17, newCigar, newPos) == 7);
300  assert(newPos == 18);
301  newCigar.getCigarString(newCigarString);
302  //std::cout << newCigarString.c_str() << std::endl;
303  assert(strcmp(newCigarString.c_str(), "8S1M") == 0);
304 
305  record.setCigar(origCigar);
306  assert(CigarHelper::softClipBeginByRefPos(record, 18, newCigar, newPos) == 8);
307  assert(newPos == 10);
308  newCigar.getCigarString(newCigarString);
309  //std::cout << newCigarString.c_str() << std::endl;
310  assert(strcmp(newCigarString.c_str(), "9S") == 0);
311 
312  record.setCigar(origCigar);
313  assert(CigarHelper::softClipBeginByRefPos(record, 19, newCigar, newPos) == 8);
314  assert(newPos == 10);
315  newCigar.getCigarString(newCigarString);
316  //std::cout << newCigarString.c_str() << std::endl;
317  assert(strcmp(newCigarString.c_str(), "9S") == 0);
318 
319  ////////////////////////////////////////////////////////
320  ////////////////////////////////////////////////////////
321  // Test clipping at every position when first non-clip instruction is delete.
322  origCigar = "3H3S3D3M3S3H";
323  record.setCigar(origCigar);
324  record.setSequence("gggAAAggg");
325  // Cigar: HHHSSSDDDMMMSSSHHH
326  // ReadPos: 000 000000
327  // ReadPos: 012 345678
328  // RefPos: 111111
329  // RefPos: 012345
330  record.setCigar(origCigar);
331  assert(CigarHelper::softClipBeginByRefPos(record, 9, newCigar, newPos) ==
332  CigarHelper::NO_CLIP);
333  assert(newPos == 10);
334  newCigar.getCigarString(newCigarString);
335  //std::cout << newCigarString.c_str() << std::endl;
336  assert(strcmp(newCigarString.c_str(), origCigar) == 0);
337 
338  record.setCigar(origCigar);
339  assert(CigarHelper::softClipBeginByRefPos(record, 10, newCigar, newPos) == 2);
340  assert(newPos == 13);
341  newCigar.getCigarString(newCigarString);
342  //std::cout << newCigarString.c_str() << std::endl;
343  assert(strcmp(newCigarString.c_str(), "3H3S3M3S3H") == 0);
344 
345  record.setCigar(origCigar);
346  assert(CigarHelper::softClipBeginByRefPos(record, 11, newCigar, newPos) == 2);
347  assert(newPos == 13);
348  newCigar.getCigarString(newCigarString);
349  //std::cout << newCigarString.c_str() << std::endl;
350  assert(strcmp(newCigarString.c_str(), "3H3S3M3S3H") == 0);
351 
352  record.setCigar(origCigar);
353  assert(CigarHelper::softClipBeginByRefPos(record, 12, newCigar, newPos) == 2);
354  assert(newPos == 13);
355  newCigar.getCigarString(newCigarString);
356  //std::cout << newCigarString.c_str() << std::endl;
357  assert(strcmp(newCigarString.c_str(), "3H3S3M3S3H") == 0);
358 
359  record.setCigar(origCigar);
360  assert(CigarHelper::softClipBeginByRefPos(record, 13, newCigar, newPos) == 3);
361  assert(newPos == 14);
362  newCigar.getCigarString(newCigarString);
363  //std::cout << newCigarString.c_str() << std::endl;
364  assert(strcmp(newCigarString.c_str(), "3H4S2M3S3H") == 0);
365 
366  record.setCigar(origCigar);
367  assert(CigarHelper::softClipBeginByRefPos(record, 14, newCigar, newPos) == 4);
368  assert(newPos == 15);
369  newCigar.getCigarString(newCigarString);
370  //std::cout << newCigarString.c_str() << std::endl;
371  assert(strcmp(newCigarString.c_str(), "3H5S1M3S3H") == 0);
372 
373  record.setCigar(origCigar);
374  assert(CigarHelper::softClipBeginByRefPos(record, 15, newCigar, newPos) == 8);
375  assert(newPos == 10);
376  newCigar.getCigarString(newCigarString);
377  //std::cout << newCigarString.c_str() << std::endl;
378  assert(strcmp(newCigarString.c_str(), "3H9S3H") == 0);
379 
380  record.setCigar(origCigar);
381  assert(CigarHelper::softClipBeginByRefPos(record, 16, newCigar, newPos) == 8);
382  assert(newPos == 10);
383  newCigar.getCigarString(newCigarString);
384  //std::cout << newCigarString.c_str() << std::endl;
385  assert(strcmp(newCigarString.c_str(), "3H9S3H") == 0);
386 
387  ////////////////////////////////////////////////////////
388  ////////////////////////////////////////////////////////
389  // Test clipping at every position when first non-clip instruction is insert.
390  origCigar = "3H3S3I3M3S3H";
391  record.setCigar(origCigar);
392  record.setSequence("gggAAATTTggg");
393  // Cigar: HHHSSSIIIMMMSSSHHH
394  // ReadPos: 000000000011
395  // ReadPos: 012345678901
396  // RefPos: 111
397  // RefPos: 012
398  record.setCigar(origCigar);
399  assert(CigarHelper::softClipBeginByRefPos(record, 9, newCigar, newPos) ==
400  CigarHelper::NO_CLIP);
401  assert(newPos == 10);
402  newCigar.getCigarString(newCigarString);
403  //std::cout << newCigarString.c_str() << std::endl;
404  assert(strcmp(newCigarString.c_str(), origCigar) == 0);
405 
406  record.setCigar(origCigar);
407  assert(CigarHelper::softClipBeginByRefPos(record, 10, newCigar, newPos) == 6);
408  assert(newPos == 11);
409  newCigar.getCigarString(newCigarString);
410  //std::cout << newCigarString.c_str() << std::endl;
411  assert(strcmp(newCigarString.c_str(), "3H7S2M3S3H") == 0);
412 
413  record.setCigar(origCigar);
414  assert(CigarHelper::softClipBeginByRefPos(record, 11, newCigar, newPos) == 7);
415  assert(newPos == 12);
416  newCigar.getCigarString(newCigarString);
417  //std::cout << newCigarString.c_str() << std::endl;
418  assert(strcmp(newCigarString.c_str(), "3H8S1M3S3H") == 0);
419 
420  record.setCigar(origCigar);
421  assert(CigarHelper::softClipBeginByRefPos(record, 12, newCigar, newPos) == 11);
422  assert(newPos == 10);
423  newCigar.getCigarString(newCigarString);
424  //std::cout << newCigarString.c_str() << std::endl;
425  assert(strcmp(newCigarString.c_str(), "3H12S3H") == 0);
426 
427  record.setCigar(origCigar);
428  assert(CigarHelper::softClipBeginByRefPos(record, 13, newCigar, newPos) == 11);
429  assert(newPos == 10);
430  newCigar.getCigarString(newCigarString);
431  //std::cout << newCigarString.c_str() << std::endl;
432  assert(strcmp(newCigarString.c_str(), "3H12S3H") == 0);
433 }
434 
435 
436 void CigarHelperTest::testSoftClipEndByRefPos()
437 {
438  SamRecord record;
439  CigarRoller newCigar;
440  std::string newCigarString;
441 
442  // Setup the current Cigar.
443  // Cigar: HHHSSSMMMDDDMMMIIIMMMPPPMMMDDDMMMSSSHHH
444  // ReadPos: 000000 000011111 111 112222
445  // ReadPos: 012345 678901234 567 890123
446  // RefPos: 111111111 122 222222223
447  // RefPos: 012345678 901 234567890
448  const char* origCigar = "3H3S3M3D3M3I3M3P3M3D3M3S3H";
449  record.setCigar(origCigar);
450  record.set0BasedPosition(10);
451  record.setSequence("gggAAATTTCCCTTTGGGAAAggg");
452 
453  ////////////////////////////////////////////////////////
454  // Clip outside of the range (after). Nothing should change.
455  assert(CigarHelper::softClipEndByRefPos(record, 10000, newCigar) ==
456  CigarHelper::NO_CLIP);
457  newCigar.getCigarString(newCigarString);
458  assert(strcmp(newCigarString.c_str(), origCigar) == 0);
459 
460  ////////////////////////////////////////////////////////
461  // Clip outside of the range (before). Everything should be clipped.
462  assert(CigarHelper::softClipEndByRefPos(record, 1, newCigar) == 0);
463  newCigar.getCigarString(newCigarString);
464  //std::cout << newCigarString.c_str() << std::endl;
465  assert(strcmp(newCigarString.c_str(), "3H24S3H") == 0);
466 
467 
468  ////////////////////////////////////////////////////////
469  ////////////////////////////////////////////////////////
470  // Test clipping at every position of the read.
471 
472  ////////////////////////////////////////////////////////
473  // Clip at the first position.
474  assert(CigarHelper::softClipEndByRefPos(record, 10, newCigar) == 0);
475  newCigar.getCigarString(newCigarString);
476  //std::cout << newCigarString.c_str() << std::endl;
477  assert(strcmp(newCigarString.c_str(), "3H24S3H") == 0);
478 
479  ////////////////////////////////////////////////////////
480  // Clip in the middle of the first Match.
481  assert(CigarHelper::softClipEndByRefPos(record, 11, newCigar) == 4);
482  newCigar.getCigarString(newCigarString);
483  //std::cout << newCigarString.c_str() << std::endl;
484  assert(strcmp(newCigarString.c_str(), "3H3S1M20S3H") == 0);
485 
486  ////////////////////////////////////////////////////////
487  // Clip just before the first deletion.
488  assert(CigarHelper::softClipEndByRefPos(record, 12, newCigar) == 5);
489  newCigar.getCigarString(newCigarString);
490  //std::cout << newCigarString.c_str() << std::endl;
491  assert(strcmp(newCigarString.c_str(), "3H3S2M19S3H") == 0);
492 
493  ////////////////////////////////////////////////////////
494  // Clip at the first deletion.
495  assert(CigarHelper::softClipEndByRefPos(record, 13, newCigar) == 6);
496  newCigar.getCigarString(newCigarString);
497  //std::cout << newCigarString.c_str() << std::endl;
498  assert(strcmp(newCigarString.c_str(), "3H3S3M18S3H") == 0);
499 
500  ////////////////////////////////////////////////////////
501  // Clip in the middle of the first deletion.
502  assert(CigarHelper::softClipEndByRefPos(record, 14, newCigar) == 6);
503  newCigar.getCigarString(newCigarString);
504  //std::cout << newCigarString.c_str() << std::endl;
505  assert(strcmp(newCigarString.c_str(), "3H3S3M18S3H") == 0);
506 
507  ////////////////////////////////////////////////////////
508  // Clip in the end of the first deletion.
509  assert(CigarHelper::softClipEndByRefPos(record, 15, newCigar) == 6);
510  newCigar.getCigarString(newCigarString);
511  //std::cout << newCigarString.c_str() << std::endl;
512  assert(strcmp(newCigarString.c_str(), "3H3S3M18S3H") == 0);
513 
514  ////////////////////////////////////////////////////////
515  // Clip just after the first deletion (should remove the deletion).
516  assert(CigarHelper::softClipEndByRefPos(record, 16, newCigar) == 6);
517  newCigar.getCigarString(newCigarString);
518  //std::cout << newCigarString.c_str() << std::endl;
519  assert(strcmp(newCigarString.c_str(), "3H3S3M18S3H") == 0);
520 
521  ////////////////////////////////////////////////////////
522  // Clip in middle of read after 1st deletion.
523  assert(CigarHelper::softClipEndByRefPos(record, 17, newCigar) == 7);
524  newCigar.getCigarString(newCigarString);
525  //std::cout << newCigarString.c_str() << std::endl;
526  assert(strcmp(newCigarString.c_str(), "3H3S3M3D1M17S3H") == 0);
527 
528  ////////////////////////////////////////////////////////
529  // Clip in middle of read after 1st deletion.
530  assert(CigarHelper::softClipEndByRefPos(record, 18, newCigar) == 8);
531  newCigar.getCigarString(newCigarString);
532  //std::cout << newCigarString.c_str() << std::endl;
533  assert(strcmp(newCigarString.c_str(), "3H3S3M3D2M16S3H") == 0);
534 
535  ////////////////////////////////////////////////////////
536  // Clip just after the 1st insertion.
537  assert(CigarHelper::softClipEndByRefPos(record, 19, newCigar) == 12);
538  newCigar.getCigarString(newCigarString);
539  //std::cout << newCigarString.c_str() << std::endl;
540  assert(strcmp(newCigarString.c_str(), "3H3S3M3D3M3I12S3H") == 0);
541 
542  ////////////////////////////////////////////////////////
543  // Clip in middle of the match after 1st insertion.
544  assert(CigarHelper::softClipEndByRefPos(record, 20, newCigar) == 13);
545  newCigar.getCigarString(newCigarString);
546  //std::cout << newCigarString.c_str() << std::endl;
547  assert(strcmp(newCigarString.c_str(), "3H3S3M3D3M3I1M11S3H") == 0);
548 
549  ////////////////////////////////////////////////////////
550  // Clip in middle of the match after 1st insertion.
551  assert(CigarHelper::softClipEndByRefPos(record, 21, newCigar) == 14);
552  newCigar.getCigarString(newCigarString);
553  //std::cout << newCigarString.c_str() << std::endl;
554  assert(strcmp(newCigarString.c_str(), "3H3S3M3D3M3I2M10S3H") == 0);
555 
556  ////////////////////////////////////////////////////////
557  // Clip right after the pad
558  assert(CigarHelper::softClipEndByRefPos(record, 22, newCigar) == 15);
559  newCigar.getCigarString(newCigarString);
560  //std::cout << newCigarString.c_str() << std::endl;
561  assert(strcmp(newCigarString.c_str(), "3H3S3M3D3M3I3M9S3H") == 0);
562 
563  ////////////////////////////////////////////////////////
564  // Clip middle of read after the pad
565  assert(CigarHelper::softClipEndByRefPos(record, 23, newCigar) == 16);
566  newCigar.getCigarString(newCigarString);
567  //std::cout << newCigarString.c_str() << std::endl;
568  assert(strcmp(newCigarString.c_str(), "3H3S3M3D3M3I3M3P1M8S3H") == 0);
569 
570  ////////////////////////////////////////////////////////
571  // Clip end of read after the pad before deletion
572  assert(CigarHelper::softClipEndByRefPos(record, 24, newCigar) == 17);
573  newCigar.getCigarString(newCigarString);
574  //std::cout << newCigarString.c_str() << std::endl;
575  assert(strcmp(newCigarString.c_str(), "3H3S3M3D3M3I3M3P2M7S3H") == 0);
576 
577  ////////////////////////////////////////////////////////
578  // Clip at start of 2nd deletion.
579  assert(CigarHelper::softClipEndByRefPos(record, 25, newCigar) == 18);
580  newCigar.getCigarString(newCigarString);
581  //std::cout << newCigarString.c_str() << std::endl;
582  assert(strcmp(newCigarString.c_str(), "3H3S3M3D3M3I3M3P3M6S3H") == 0);
583 
584  ////////////////////////////////////////////////////////
585  // Clip in 2nd deletion.
586  assert(CigarHelper::softClipEndByRefPos(record, 26, newCigar) == 18);
587  newCigar.getCigarString(newCigarString);
588  //std::cout << newCigarString.c_str() << std::endl;
589  assert(strcmp(newCigarString.c_str(), "3H3S3M3D3M3I3M3P3M6S3H") == 0);
590 
591  ////////////////////////////////////////////////////////
592  // Clip in 2nd deletion.
593  assert(CigarHelper::softClipEndByRefPos(record, 27, newCigar) == 18);
594  newCigar.getCigarString(newCigarString);
595  //std::cout << newCigarString.c_str() << std::endl;
596  assert(strcmp(newCigarString.c_str(), "3H3S3M3D3M3I3M3P3M6S3H") == 0);
597 
598  ////////////////////////////////////////////////////////
599  // Clip right after 2nd deletion.
600  assert(CigarHelper::softClipEndByRefPos(record, 28, newCigar) == 18);
601  newCigar.getCigarString(newCigarString);
602  //std::cout << newCigarString.c_str() << std::endl;
603  assert(strcmp(newCigarString.c_str(), "3H3S3M3D3M3I3M3P3M6S3H") == 0);
604 
605  ////////////////////////////////////////////////////////
606  // Clip in middle of last match.
607  assert(CigarHelper::softClipEndByRefPos(record, 29, newCigar) == 19);
608  newCigar.getCigarString(newCigarString);
609  //std::cout << newCigarString.c_str() << std::endl;
610  assert(strcmp(newCigarString.c_str(), "3H3S3M3D3M3I3M3P3M3D1M5S3H") == 0);
611 
612  ////////////////////////////////////////////////////////
613  // Clip in middle of last match.
614  assert(CigarHelper::softClipEndByRefPos(record, 30, newCigar) == 20);
615  newCigar.getCigarString(newCigarString);
616  //std::cout << newCigarString.c_str() << std::endl;
617  assert(strcmp(newCigarString.c_str(), "3H3S3M3D3M3I3M3P3M3D2M4S3H") == 0);
618 
619  ////////////////////////////////////////////////////////
620  // Clip right after the read (no change).
621  assert(CigarHelper::softClipEndByRefPos(record, 31, newCigar) ==
622  CigarHelper::NO_CLIP);
623  newCigar.getCigarString(newCigarString);
624  //std::cout << newCigarString.c_str() << std::endl;
625  assert(strcmp(newCigarString.c_str(), origCigar) == 0);
626 
627  ////////////////////////////////////////////////////////
628  ////////////////////////////////////////////////////////
629  // Test clipping at every position when insertions & deletions
630  // are next to each other.
631  origCigar = "3M3D3I3M";
632  record.setCigar(origCigar);
633  record.setSequence("GGGAAAGGG");
634  // Cigar: MMMDDDIIIMMM
635  // ReadPos: 000 000000
636  // ReadPos: 012 345678
637  // RefPos: 111111 111
638  // RefPos: 012345 678
639  record.setCigar(origCigar);
640  assert(CigarHelper::softClipEndByRefPos(record, 9, newCigar) == 0);
641  newCigar.getCigarString(newCigarString);
642  //std::cout << newCigarString.c_str() << std::endl;
643  assert(strcmp(newCigarString.c_str(), "9S") == 0);
644 
645  record.setCigar(origCigar);
646  assert(CigarHelper::softClipEndByRefPos(record, 10, newCigar) == 0);
647  newCigar.getCigarString(newCigarString);
648  //std::cout << newCigarString.c_str() << std::endl;
649  assert(strcmp(newCigarString.c_str(), "9S") == 0);
650 
651  record.setCigar(origCigar);
652  assert(CigarHelper::softClipEndByRefPos(record, 11, newCigar) == 1);
653  newCigar.getCigarString(newCigarString);
654  //std::cout << newCigarString.c_str() << std::endl;
655  assert(strcmp(newCigarString.c_str(), "1M8S") == 0);
656 
657  record.setCigar(origCigar);
658  assert(CigarHelper::softClipEndByRefPos(record, 12, newCigar) == 2);
659  newCigar.getCigarString(newCigarString);
660  //std::cout << newCigarString.c_str() << std::endl;
661  assert(strcmp(newCigarString.c_str(), "2M7S") == 0);
662 
663  record.setCigar(origCigar);
664  assert(CigarHelper::softClipEndByRefPos(record, 13, newCigar) == 3);
665  newCigar.getCigarString(newCigarString);
666  //std::cout << newCigarString.c_str() << std::endl;
667  assert(strcmp(newCigarString.c_str(), "3M6S") == 0);
668 
669  record.setCigar(origCigar);
670  assert(CigarHelper::softClipEndByRefPos(record, 14, newCigar) == 3);
671  newCigar.getCigarString(newCigarString);
672  //std::cout << newCigarString.c_str() << std::endl;
673  assert(strcmp(newCigarString.c_str(), "3M6S") == 0);
674 
675  record.setCigar(origCigar);
676  assert(CigarHelper::softClipEndByRefPos(record, 15, newCigar) == 3);
677  newCigar.getCigarString(newCigarString);
678  //std::cout << newCigarString.c_str() << std::endl;
679  assert(strcmp(newCigarString.c_str(), "3M6S") == 0);
680 
681  record.setCigar(origCigar);
682  assert(CigarHelper::softClipEndByRefPos(record, 16, newCigar) == 6);
683  newCigar.getCigarString(newCigarString);
684  //std::cout << newCigarString.c_str() << std::endl;
685  assert(strcmp(newCigarString.c_str(), "3M3D3I3S") == 0);
686 
687  record.setCigar(origCigar);
688  assert(CigarHelper::softClipEndByRefPos(record, 17, newCigar) == 7);
689  newCigar.getCigarString(newCigarString);
690  //std::cout << newCigarString.c_str() << std::endl;
691  assert(strcmp(newCigarString.c_str(), "3M3D3I1M2S") == 0);
692 
693  record.setCigar(origCigar);
694  assert(CigarHelper::softClipEndByRefPos(record, 18, newCigar) == 8);
695  newCigar.getCigarString(newCigarString);
696  //std::cout << newCigarString.c_str() << std::endl;
697  assert(strcmp(newCigarString.c_str(), "3M3D3I2M1S") == 0);
698 
699  record.setCigar(origCigar);
700  assert(CigarHelper::softClipEndByRefPos(record, 19, newCigar) ==
701  CigarHelper::NO_CLIP);
702  newCigar.getCigarString(newCigarString);
703  //std::cout << newCigarString.c_str() << std::endl;
704  assert(strcmp(newCigarString.c_str(), origCigar) == 0);
705 
706  ////////////////////////////////////////////////////////
707  ////////////////////////////////////////////////////////
708  // Test clipping at every position when first non-clip instruction is delete.
709  origCigar = "3H3S3D3M3S3H";
710  record.setCigar(origCigar);
711  record.setSequence("gggAAAggg");
712  // Cigar: HHHSSSDDDMMMSSSHHH
713  // ReadPos: 000 000000
714  // ReadPos: 012 345678
715  // RefPos: 111111
716  // RefPos: 012345
717  record.setCigar(origCigar);
718  assert(CigarHelper::softClipEndByRefPos(record, 9, newCigar) == 0);
719  newCigar.getCigarString(newCigarString);
720  //std::cout << newCigarString.c_str() << std::endl;
721  assert(strcmp(newCigarString.c_str(), "3H9S3H") == 0);
722 
723  record.setCigar(origCigar);
724  assert(CigarHelper::softClipEndByRefPos(record, 10, newCigar) == 0);
725  newCigar.getCigarString(newCigarString);
726  //std::cout << newCigarString.c_str() << std::endl;
727  assert(strcmp(newCigarString.c_str(), "3H9S3H") == 0);
728 
729  record.setCigar(origCigar);
730  assert(CigarHelper::softClipEndByRefPos(record, 11, newCigar) == 0);
731  newCigar.getCigarString(newCigarString);
732  //std::cout << newCigarString.c_str() << std::endl;
733  assert(strcmp(newCigarString.c_str(), "3H9S3H") == 0);
734 
735  record.setCigar(origCigar);
736  assert(CigarHelper::softClipEndByRefPos(record, 12, newCigar) == 0);
737  newCigar.getCigarString(newCigarString);
738  //std::cout << newCigarString.c_str() << std::endl;
739  assert(strcmp(newCigarString.c_str(), "3H9S3H") == 0);
740 
741  record.setCigar(origCigar);
742  assert(CigarHelper::softClipEndByRefPos(record, 13, newCigar) == 0);
743  newCigar.getCigarString(newCigarString);
744  //std::cout << newCigarString.c_str() << std::endl;
745  assert(strcmp(newCigarString.c_str(), "3H9S3H") == 0);
746 
747  record.setCigar(origCigar);
748  assert(CigarHelper::softClipEndByRefPos(record, 14, newCigar) == 4);
749  newCigar.getCigarString(newCigarString);
750  //std::cout << newCigarString.c_str() << std::endl;
751  assert(strcmp(newCigarString.c_str(), "3H3S3D1M5S3H") == 0);
752 
753  record.setCigar(origCigar);
754  assert(CigarHelper::softClipEndByRefPos(record, 15, newCigar) == 5);
755  newCigar.getCigarString(newCigarString);
756  //std::cout << newCigarString.c_str() << std::endl;
757  assert(strcmp(newCigarString.c_str(), "3H3S3D2M4S3H") == 0);
758 
759  record.setCigar(origCigar);
760  assert(CigarHelper::softClipEndByRefPos(record, 16, newCigar) ==
761  CigarHelper::NO_CLIP);
762  newCigar.getCigarString(newCigarString);
763  //std::cout << newCigarString.c_str() << std::endl;
764  assert(strcmp(newCigarString.c_str(), origCigar) == 0);
765 
766  ////////////////////////////////////////////////////////
767  ////////////////////////////////////////////////////////
768  // Test clipping at every position when first non-clip instruction is insert.
769  origCigar = "3H3S3I3M3S3H";
770  record.setCigar(origCigar);
771  record.setSequence("gggAAATTTggg");
772  // Cigar: HHHSSSIIIMMMSSSHHH
773  // ReadPos: 000000000011
774  // ReadPos: 012345678901
775  // RefPos: 111
776  // RefPos: 012
777  record.setCigar(origCigar);
778  assert(CigarHelper::softClipEndByRefPos(record, 9, newCigar) == 0);
779  newCigar.getCigarString(newCigarString);
780  //std::cout << newCigarString.c_str() << std::endl;
781  assert(strcmp(newCigarString.c_str(), "3H12S3H") == 0);
782 
783  record.setCigar(origCigar);
784  assert(CigarHelper::softClipEndByRefPos(record, 10, newCigar) == 6);
785  newCigar.getCigarString(newCigarString);
786  //std::cout << newCigarString.c_str() << std::endl;
787  assert(strcmp(newCigarString.c_str(), "3H3S3I6S3H") == 0);
788 
789  record.setCigar(origCigar);
790  assert(CigarHelper::softClipEndByRefPos(record, 11, newCigar) == 7);
791  newCigar.getCigarString(newCigarString);
792  //std::cout << newCigarString.c_str() << std::endl;
793  assert(strcmp(newCigarString.c_str(), "3H3S3I1M5S3H") == 0);
794 
795  record.setCigar(origCigar);
796  assert(CigarHelper::softClipEndByRefPos(record, 12, newCigar) == 8);
797  newCigar.getCigarString(newCigarString);
798  //std::cout << newCigarString.c_str() << std::endl;
799  assert(strcmp(newCigarString.c_str(), "3H3S3I2M4S3H") == 0);
800 
801  record.setCigar(origCigar);
802  assert(CigarHelper::softClipEndByRefPos(record, 13, newCigar) ==
803  CigarHelper::NO_CLIP);
804  newCigar.getCigarString(newCigarString);
805  //std::cout << newCigarString.c_str() << std::endl;
806  assert(strcmp(newCigarString.c_str(), origCigar) == 0);
807 }
static int32_t softClipBeginByRefPos(SamRecord &record, int32_t refPosition0Based, CigarRoller &newCigar, int32_t &new0BasedPosition)
Soft clip the cigar from the beginning of the read at the specified reference position.
Definition: CigarHelper.cpp:23
static int32_t softClipEndByRefPos(SamRecord &record, int32_t refPosition0Based, CigarRoller &newCigar)
Soft clip the cigar from the back of the read at the specified reference position.
The purpose of this class is to provide accessors for setting, updating, modifying the CIGAR object....
Definition: CigarRoller.h:67
void getCigarString(String &cigarString) const
Set the passed in String to the string reprentation of the Cigar operations in this object.
Definition: Cigar.cpp:52
Class providing an easy to use interface to get/set/operate on the fields in a SAM/BAM record.
Definition: SamRecord.h:52
bool setCigar(const char *cigar)
Set the CIGAR to the specified SAM formatted cigar string.
Definition: SamRecord.cpp:259
bool setSequence(const char *seq)
Sets the sequence (SEQ) to the specified SAM formatted sequence string.
Definition: SamRecord.cpp:344
bool set0BasedPosition(int32_t position)
Set the leftmost position using the specified 0-based (BAM format) value.
Definition: SamRecord.cpp:242