libStatGen Software 1
EqualsTest Class Reference

Public Types

enum  FileType { SAM , BAM }
 

Static Public Member Functions

static void testEq (FileType inputType)
 

Detailed Description

Definition at line 21 of file TestEquals.h.

Member Enumeration Documentation

◆ FileType

enum EqualsTest::FileType

Definition at line 24 of file TestEquals.h.

24{SAM, BAM};

Member Function Documentation

◆ testEq()

void EqualsTest::testEq ( FileType  inputType)
static

Definition at line 71 of file TestEquals.cpp.

72{
73 reset();
74 SamFile inSam;
75
76 std::string outputBase = "results/out";
77
78 if(inputType == SAM)
79 {
80 assert(inSam.OpenForRead("testFiles/testEq.sam"));
81 outputBase += "SamEq";
82 }
83 else
84 {
85 assert(inSam.OpenForRead("testFiles/testEq.bam"));
86 outputBase += "BamEq";
87 }
88
89 // Read the SAM Header.
90 SamFileHeader samHeader;
91 assert(inSam.ReadHeader(samHeader));
92
93 std::string outputName = outputBase + "Bases.sam";
94 SamFile outBasesSam( outputName.c_str(), SamFile::WRITE);
95 outputName = outputBase + "Equals.sam";
96 SamFile outEqualsSam(outputName.c_str(), SamFile::WRITE);
97 outputName = outputBase + "Orig.sam";
98 SamFile outOrigSam( outputName.c_str(), SamFile::WRITE);
99 outputName = outputBase + "Bases.bam";
100 SamFile outBasesBam( outputName.c_str(), SamFile::WRITE);
101 outputName = outputBase + "Equals.bam";
102 SamFile outEqualsBam(outputName.c_str(), SamFile::WRITE);
103 outputName = outputBase + "Orig.bam";
104 SamFile outOrigBam( outputName.c_str(), SamFile::WRITE);
105 assert(outBasesSam.WriteHeader(samHeader));
106 assert(outEqualsSam.WriteHeader(samHeader));
107 assert(outOrigSam.WriteHeader(samHeader));
108 assert(outBasesBam.WriteHeader(samHeader));
109 assert(outEqualsBam.WriteHeader(samHeader));
110 assert(outOrigBam.WriteHeader(samHeader));
111
112 outBasesSam.SetWriteSequenceTranslation(SamRecord::BASES);
113 outEqualsSam.SetWriteSequenceTranslation(SamRecord::EQUAL);
114 outOrigSam.SetWriteSequenceTranslation(SamRecord::NONE);
115 outBasesBam.SetWriteSequenceTranslation(SamRecord::BASES);
116 outEqualsBam.SetWriteSequenceTranslation(SamRecord::EQUAL);
117 outOrigBam.SetWriteSequenceTranslation(SamRecord::NONE);
118
119 GenomeSequence reference("testFiles/chr1_partial.fa");
120
121 inSam.SetReference(&reference);
122
123 SamRecord samRecord;
124
125 // The set of 16 variations are repeated 3 times: once with all charcters
126 // 1) Matches have the actual bases in them.
127 // 2) Matches have '='
128 // 3) Matches are mixed between bases and '='
129 // Since Sequences are 4 characters long, there are 16 variations
130 // of match/mismatch.
131 for(int j = 0; j < 16; j++)
132 {
133 assert(inSam.ReadRecord(samHeader, samRecord) == true);
134 validateEqRead(samRecord, j, READ_SEQS_BASES[j]);
135 assert(outBasesSam.WriteRecord(samHeader, samRecord));
136 assert(outEqualsSam.WriteRecord(samHeader, samRecord));
137 assert(outOrigSam.WriteRecord(samHeader, samRecord));
138 assert(outBasesBam.WriteRecord(samHeader, samRecord));
139 assert(outEqualsBam.WriteRecord(samHeader, samRecord));
140 assert(outOrigBam.WriteRecord(samHeader, samRecord));
141 }
142 for(int j = 0; j < 16; j++)
143 {
144 assert(inSam.ReadRecord(samHeader, samRecord) == true);
145 validateEqRead(samRecord, j, READ_SEQS_EQUALS[j]);
146 assert(outBasesSam.WriteRecord(samHeader, samRecord));
147 assert(outEqualsSam.WriteRecord(samHeader, samRecord));
148 assert(outOrigSam.WriteRecord(samHeader, samRecord));
149 assert(outBasesBam.WriteRecord(samHeader, samRecord));
150 assert(outEqualsBam.WriteRecord(samHeader, samRecord));
151 assert(outOrigBam.WriteRecord(samHeader, samRecord));
152 }
153 for(int j = 0; j < 16; j++)
154 {
155 assert(inSam.ReadRecord(samHeader, samRecord) == true);
156 validateEqRead(samRecord, j, READ_SEQS_MIXED[j]);
157 assert(outBasesSam.WriteRecord(samHeader, samRecord));
158 assert(outEqualsSam.WriteRecord(samHeader, samRecord));
159 assert(outOrigSam.WriteRecord(samHeader, samRecord));
160 assert(outBasesBam.WriteRecord(samHeader, samRecord));
161 assert(outEqualsBam.WriteRecord(samHeader, samRecord));
162 assert(outOrigBam.WriteRecord(samHeader, samRecord));
163 }
164
165 expectedCigar = "2S4M";
166 expectedCigarHex.clear();
167 expectedCigarHex.push_back(0x24);
168 expectedCigarHex.push_back(0x40);
169 expected0BasedUnclippedStart = expectedRecord.myPosition-2;
170 expected1BasedUnclippedStart = expected0BasedUnclippedStart + 1;
171 expectedRecord.myBlockSize = 70;
172 expectedRecord.myReadNameLength = 21;
173 expectedRecord.myCigarLength = 2;
174 expectedRecord.myReadLength = 6;
175 expectedQuality = "??I00?";
176 assert(inSam.ReadRecord(samHeader, samRecord) == true);
177 validateEqRead(samRecord, 16, READ_SEQS_MIXED[16]);
178 assert(outBasesSam.WriteRecord(samHeader, samRecord));
179 assert(outEqualsSam.WriteRecord(samHeader, samRecord));
180 assert(outOrigSam.WriteRecord(samHeader, samRecord));
181 assert(outBasesBam.WriteRecord(samHeader, samRecord));
182 assert(outEqualsBam.WriteRecord(samHeader, samRecord));
183 assert(outOrigBam.WriteRecord(samHeader, samRecord));
184
185 expectedCigar = "4M4H";
186 expectedCigarHex.clear();
187 expectedCigarHex.push_back(0x40);
188 expectedCigarHex.push_back(0x45);
189 expected0BasedUnclippedStart = expectedRecord.myPosition;
190 expected1BasedUnclippedStart = expected0BasedUnclippedStart + 1;
191 expected0BasedUnclippedEnd = expectedRecord.myPosition + 7;
192 expected1BasedUnclippedEnd = expected0BasedUnclippedEnd + 1;
193 expectedRecord.myBlockSize = 65;
194 expectedRecord.myReadNameLength = 19;
195 expectedRecord.myCigarLength = 2;
196 expectedRecord.myReadLength = 4;
197 expectedQuality = "I00?";
198 assert(inSam.ReadRecord(samHeader, samRecord) == true);
199 validateEqRead(samRecord, 17, READ_SEQS_MIXED[17]);
200 assert(outBasesSam.WriteRecord(samHeader, samRecord));
201 assert(outEqualsSam.WriteRecord(samHeader, samRecord));
202 assert(outOrigSam.WriteRecord(samHeader, samRecord));
203 assert(outBasesBam.WriteRecord(samHeader, samRecord));
204 assert(outEqualsBam.WriteRecord(samHeader, samRecord));
205 assert(outOrigBam.WriteRecord(samHeader, samRecord));
206
207 expectedCigar = "1M1P1M1I1M3D1M";
208 expectedCigarHex.clear();
209 expectedCigarHex.push_back(0x10);
210 expectedCigarHex.push_back(0x16);
211 expectedCigarHex.push_back(0x10);
212 expectedCigarHex.push_back(0x11);
213 expectedCigarHex.push_back(0x10);
214 expectedCigarHex.push_back(0x32);
215 expectedCigarHex.push_back(0x10);
216 expected0BasedAlignmentEnd = expectedRecord.myPosition + 6;
217 expected1BasedAlignmentEnd = expected0BasedAlignmentEnd + 1;
218 expectedAlignmentLength = 7;
219 expected0BasedUnclippedStart = expectedRecord.myPosition;
220 expected1BasedUnclippedStart = expected0BasedUnclippedStart + 1;
221 expected0BasedUnclippedEnd = expected0BasedAlignmentEnd;
222 expected1BasedUnclippedEnd = expected0BasedUnclippedEnd + 1;
223 expectedRecord.myBlockSize = 95;
224 expectedRecord.myReadNameLength = 27;
225 expectedRecord.myCigarLength = 7;
226 expectedRecord.myReadLength = 5;
227 expectedQuality = "I00??";
228 assert(inSam.ReadRecord(samHeader, samRecord) == true);
229 validateEqRead(samRecord, 18, READ_SEQS_MIXED[18]);
230 assert(outBasesSam.WriteRecord(samHeader, samRecord));
231 assert(outEqualsSam.WriteRecord(samHeader, samRecord));
232 assert(outOrigSam.WriteRecord(samHeader, samRecord));
233 assert(outBasesBam.WriteRecord(samHeader, samRecord));
234 assert(outEqualsBam.WriteRecord(samHeader, samRecord));
235 assert(outOrigBam.WriteRecord(samHeader, samRecord));
236
237 expectedCigar = "2M2N2M";
238 expectedCigarHex.clear();
239 expectedCigarHex.push_back(0x20);
240 expectedCigarHex.push_back(0x23);
241 expectedCigarHex.push_back(0x20);
242 expected0BasedAlignmentEnd = expectedRecord.myPosition + 5;
243 expected1BasedAlignmentEnd = expected0BasedAlignmentEnd + 1;
244 expectedAlignmentLength = 6;
245 expected0BasedUnclippedStart = expectedRecord.myPosition;
246 expected1BasedUnclippedStart = expected0BasedUnclippedStart + 1;
247 expected0BasedUnclippedEnd = expected0BasedAlignmentEnd;
248 expected1BasedUnclippedEnd = expected0BasedUnclippedEnd + 1;
249 expectedRecord.myBlockSize = 74;
250 expectedRecord.myReadNameLength = 24;
251 expectedRecord.myCigarLength = 3;
252 expectedRecord.myReadLength = 4;
253 expectedQuality = "I00?";
254 assert(inSam.ReadRecord(samHeader, samRecord) == true);
255 validateEqRead(samRecord, 19, READ_SEQS_MIXED[19]);
256 assert(outBasesSam.WriteRecord(samHeader, samRecord));
257 assert(outEqualsSam.WriteRecord(samHeader, samRecord));
258 assert(outOrigSam.WriteRecord(samHeader, samRecord));
259 assert(outBasesBam.WriteRecord(samHeader, samRecord));
260 assert(outEqualsBam.WriteRecord(samHeader, samRecord));
261 assert(outOrigBam.WriteRecord(samHeader, samRecord));
262
263 // Test getNextMatchMismatch.
264 SamSingleBaseMatchInfo matchTest;
265 SamQuerySeqWithRefIter queryIter(samRecord, reference, true);
266 assert(queryIter.getNextMatchMismatch(matchTest) == true);
267 assert(matchTest.getType() == SamSingleBaseMatchInfo::MATCH);
268 assert(matchTest.getQueryIndex() == 0);
269 assert(queryIter.getNextMatchMismatch(matchTest) == true);
270 assert(matchTest.getType() == SamSingleBaseMatchInfo::MATCH);
271 assert(matchTest.getQueryIndex() == 1);
272 assert(queryIter.getNextMatchMismatch(matchTest) == true);
273 assert(matchTest.getType() == SamSingleBaseMatchInfo::MATCH);
274 assert(matchTest.getQueryIndex() == 2);
275 assert(queryIter.getNextMatchMismatch(matchTest) == true);
276 assert(matchTest.getType() == SamSingleBaseMatchInfo::MATCH);
277 assert(matchTest.getQueryIndex() == 3);
278 assert(queryIter.getNextMatchMismatch(matchTest) == false);
279
280 // Check the read that is on a different chormosome not
281 // found in the reference.
282 reset();
283 expectedRecord.myBlockSize = 56;
284 expectedRecord.myReadNameLength = 14;
285 expectedRecord.myReferenceID = 1;
286 expectedReferenceName = "2";
287 expectedRecord.myMateReferenceID = 1;
288 expectedMateReferenceName = "2";
289
290 assert(inSam.ReadRecord(samHeader, samRecord) == true);
291 validateEqRead(samRecord, 20, READ_SEQS_MIXED[20]);
292 assert(outBasesSam.WriteRecord(samHeader, samRecord));
293 assert(outEqualsSam.WriteRecord(samHeader, samRecord));
294 assert(outOrigSam.WriteRecord(samHeader, samRecord));
295 assert(outBasesBam.WriteRecord(samHeader, samRecord));
296 assert(outEqualsBam.WriteRecord(samHeader, samRecord));
297 assert(outOrigBam.WriteRecord(samHeader, samRecord));
298
299 // Check the read that is on a different chormosome and
300 // has '=', but the chromosome is not found in the reference.
301 reset();
302 expectedRecord.myBlockSize = 57;
303 expectedRecord.myReadNameLength = 15;
304 expectedRecord.myReferenceID = 1;
305 expectedReferenceName = "2";
306 expectedRecord.myMateReferenceID = 1;
307 expectedMateReferenceName = "2";
308
309 assert(inSam.ReadRecord(samHeader, samRecord) == true);
310 validateEqRead(samRecord, 21, READ_SEQS_MIXED[21]);
311 assert(outBasesSam.WriteRecord(samHeader, samRecord));
312 assert(outEqualsSam.WriteRecord(samHeader, samRecord));
313 assert(outOrigSam.WriteRecord(samHeader, samRecord));
314 assert(outBasesBam.WriteRecord(samHeader, samRecord));
315 assert(outEqualsBam.WriteRecord(samHeader, samRecord));
316 assert(outOrigBam.WriteRecord(samHeader, samRecord));
317
318 SamQuerySeqWithRefIter queryIter2(samRecord, reference, true);
319 assert(queryIter2.getNextMatchMismatch(matchTest) == false);
320}
Create/Access/Modify/Load Genome Sequences stored as binary mapped files.
This class allows a user to get/set the fields in a SAM/BAM Header.
Definition: SamFileHeader.h:35
Allows the user to easily read/write a SAM/BAM file.
Definition: SamFile.h:36
bool ReadHeader(SamFileHeader &header)
Reads the header section from the file and stores it in the passed in header.
Definition: SamFile.cpp:450
bool ReadRecord(SamFileHeader &header, SamRecord &record)
Reads the next record from the file & stores it in the passed in record.
Definition: SamFile.cpp:514
void SetReference(GenomeSequence *reference)
Sets the reference to the specified genome sequence object.
Definition: SamFile.cpp:380
@ WRITE
open for writing.
Definition: SamFile.h:41
bool OpenForRead(const char *filename, SamFileHeader *header=NULL)
Open a sam/bam file for reading with the specified filename, determing the type of file and SAM/BAM b...
Definition: SamFile.cpp:93
Iterates through the query and compare with reference.
Class providing an easy to use interface to get/set/operate on the fields in a SAM/BAM record.
Definition: SamRecord.h:52
@ NONE
Leave the sequence as is.
Definition: SamRecord.h:58
@ BASES
Translate '=' to the actual base.
Definition: SamRecord.h:60
@ EQUAL
Translate bases that match the reference to '='.
Definition: SamRecord.h:59
This class contains the match/mismatch information between the reference and a read for a single base...
int32_t getQueryIndex()
Get the query index for this object.
Type getType()
Get the type (match/mismatch/unknown) for this object.

The documentation for this class was generated from the following files: