ReadFiles.cpp

00001 /*
00002  *  Copyright (C) 2010  Regents of the University of Michigan
00003  *
00004  *   This program is free software: you can redistribute it and/or modify
00005  *   it under the terms of the GNU General Public License as published by
00006  *   the Free Software Foundation, either version 3 of the License, or
00007  *   (at your option) any later version.
00008  *
00009  *   This program is distributed in the hope that it will be useful,
00010  *   but WITHOUT ANY WARRANTY; without even the implied warranty of
00011  *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00012  *   GNU General Public License for more details.
00013  *
00014  *   You should have received a copy of the GNU General Public License
00015  *   along with this program.  If not, see <http://www.gnu.org/licenses/>.
00016  */
00017 
00018 #include "ReadFiles.h"
00019 #include "Validate.h"
00020 #include "SamTags.h"
00021 #include <assert.h>
00022 
00023 void testReadSam()
00024 {
00025     SamFile inSam;
00026     assert(inSam.OpenForRead("testFiles/testSam.sam"));
00027 
00028     // Call generic test which since the sam and bam are identical, should
00029     // contain the same results.
00030     testRead(inSam);
00031 }
00032 
00033 void testReadBam()
00034 {
00035     SamFile inSam;
00036     assert(inSam.OpenForRead("testFiles/testBam.bam"));
00037 
00038     // Call generic test which since the sam and bam are identical, should
00039     // contain the same results.
00040     testRead(inSam);
00041 }
00042 
00043 void testRead(SamFile &inSam)
00044 {
00045     // Read the SAM Header.
00046     SamFileHeader samHeader;
00047     assert(inSam.ReadHeader(samHeader));
00048 
00049     validateHeader(samHeader);
00050     testModHeader(samHeader);
00051 
00052     SamRecord samRecord;
00053     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00054     validateRead1(samRecord);
00055 
00056     // Set a new quality and get the buffer.
00057     samRecord.setQuality("ABCDE");
00058     validateRead1ModQuality(samRecord);
00059     //   void* buffer = samRecord.getRecordBuffer();
00060 
00061     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00062     validateRead2(samRecord);
00063 
00064     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00065     validateRead3(samRecord);
00066 
00067     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00068     validateRead4(samRecord);
00069 
00070     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00071     validateRead5(samRecord);
00072 
00073     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00074     validateRead6(samRecord);
00075 
00076     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00077     validateRead7(samRecord);
00078 
00079     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00080     validateRead8(samRecord);
00081 
00082     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00083     validateRead9(samRecord);
00084 
00085     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00086     validateRead10(samRecord);
00087 }
00088 
00089 
00090 void testAddHeaderAndTagToFile(const char* inputName, const char* outputName)
00091 {
00092     SamFile inSam, outSam;
00093     assert(inSam.OpenForRead(inputName));
00094     assert(outSam.OpenForWrite(outputName));
00095 
00096     // Read the SAM Header.
00097     SamFileHeader samHeader;
00098     assert(inSam.ReadHeader(samHeader));
00099 
00100     // Add a header line.
00101     assert(samHeader.addHeaderLine("@RG\tID:myID\tSM:mySM") == false);
00102     assert(samHeader.addHeaderLine("@RG\tID:myID3\tSM:mySM") == true);
00103 
00104     // Write Header
00105     assert(outSam.WriteHeader(samHeader));
00106 
00107     SamRecord samRecord;
00108     assert(inSam.ReadRecord(samHeader, samRecord));
00109     //   validateRead1(samRecord);
00110     // Add two tags.
00111     assert(samRecord.addIntTag("XA", 123));
00112     assert(samRecord.addIntTag("XA", 456));
00113     assert(samRecord.addTag("RR", 'Z', "myID1"));
00114     assert(samRecord.addTag("RR", 'Z', "myID2"));
00115 
00116     // Write as Sam.
00117     assert(outSam.WriteRecord(samHeader, samRecord));
00118 
00119     // TODO, add test to verify it was written correctly.
00120 
00121     // Read a couple of records to make sure it properly can read them even
00122     // if they are bigger than the original.
00123     assert(inSam.ReadRecord(samHeader, samRecord));
00124     assert(inSam.ReadRecord(samHeader, samRecord));
00125 
00126     //  Check the MD tag, which requires the reference.
00127     GenomeSequence reference("testFiles/chr1_partial.fa");
00128     assert(SamTags::isMDTagCorrect(samRecord, reference) == false);
00129     String newMDTag;
00130     SamTags::createMDTag(newMDTag, samRecord, reference);
00131     assert(newMDTag == "2T1N0");
00132     assert(SamTags::updateMDTag(samRecord, reference));
00133     // Write as Sam.
00134     assert(outSam.WriteRecord(samHeader, samRecord));
00135 }
00136 
00137 
00138 // Test reading a file, validating it is sorted.
00139 void testValidateSortedRead()
00140 {
00141     // Open a file for reading.
00142     SamFile inSam(ErrorHandler::RETURN);
00143     assert(inSam.OpenForRead("testFiles/testSam.sam"));
00144 
00145     // Set the validation to COORDINATE.
00146     inSam.setSortedValidation(SamFile::COORDINATE);
00147 
00148     // Read the SAM Header.
00149     SamFileHeader samHeader;
00150     assert(inSam.ReadHeader(samHeader));
00151     
00152     SamRecord samRecord;
00153     // Succeed, first record.
00154     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00155     validateRead1(samRecord);
00156 
00157     // Succeed, higher coordinate.
00158     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00159     validateRead2(samRecord);
00160 
00161     // Failed sort order - due to coord.
00162     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00163     validateRead3(samRecord);
00164 
00165     // Failed sort order - due to coord.
00166     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00167     validateRead4(samRecord);
00168 
00169     // Succeed, new reference id
00170     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00171     validateRead5(samRecord);
00172 
00173     // Fail, previous reference id.
00174     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00175     validateRead6(samRecord);
00176 
00177     // Succeed, same reference id, higher coord.
00178     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00179     validateRead7(samRecord);
00180 
00181     // Succeed, *, new reference id.
00182     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00183     validateRead8(samRecord);
00184 
00185     // Fail, reference id is not *
00186     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00187     validateRead9(samRecord);
00188 
00189     // Succeed, valid reference id, and no coordinate.
00190     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00191     validateRead10(samRecord);
00192 
00193 
00194     ////////////////////////////////////////////
00195     // Reopen the file for reading
00196     assert(inSam.OpenForRead("testFiles/testSam.sam"));
00197 
00198     // Set the validation to QUERY_NAME.
00199     inSam.setSortedValidation(SamFile::QUERY_NAME);
00200 
00201     // Read the SAM Header.
00202     assert(inSam.ReadHeader(samHeader));
00203    
00204     // Succeed, first record.
00205     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00206     validateRead1(samRecord);
00207 
00208     // Succeed, same name.
00209     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00210     validateRead2(samRecord);
00211 
00212     // Failed sort order - due to alphabetical
00213     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00214     validateRead3(samRecord);
00215 
00216     // Failed sort order - due to alphabetical
00217     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00218     validateRead4(samRecord);
00219 
00220     // Failed sort order - due to alphabetical
00221     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00222     validateRead5(samRecord);
00223 
00224     // Failed sort order - due to alphabetical
00225     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00226     validateRead6(samRecord);
00227 
00228     // Failed sort order - due to alphabetical
00229     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00230     validateRead7(samRecord);
00231 
00232     // Succeed
00233     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00234     validateRead8(samRecord);
00235 
00236     // Failed sort order - due to alphabetical
00237     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00238     validateRead9(samRecord);
00239 
00240     // Succeed
00241     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00242     validateRead10(samRecord);
00243 
00244     ////////////////////////////////////////////
00245     // Reopen the file for reading
00246     assert(inSam.OpenForRead("testFiles/testSam.sam"));
00247 
00248     // Set the validation to the SO Flag.  Not set, so it is UNSORTED, so 
00249     // all reads should pass.
00250     inSam.setSortedValidation(SamFile::FLAG);
00251 
00252     // Read the SAM Header.
00253     assert(inSam.ReadHeader(samHeader));
00254    
00255     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00256     validateRead1(samRecord);
00257 
00258     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00259     validateRead2(samRecord);
00260 
00261     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00262     validateRead3(samRecord);
00263 
00264     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00265     validateRead4(samRecord);
00266 
00267     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00268     validateRead5(samRecord);
00269 
00270     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00271     validateRead6(samRecord);
00272 
00273     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00274     validateRead7(samRecord);
00275 
00276     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00277     validateRead8(samRecord);
00278 
00279     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00280     validateRead9(samRecord);
00281 
00282     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00283     validateRead10(samRecord);
00284 
00285     ////////////////////////////////////////////
00286     // Reopen for reading SO FLAG set to coordinate.
00287     assert(inSam.OpenForRead("testFiles/testSamSOcoord.sam"));
00288 
00289     // Set the validation to SO FLAG which is set to coordinate.
00290     inSam.setSortedValidation(SamFile::FLAG);
00291 
00292     // Read the SAM Header.
00293     assert(inSam.ReadHeader(samHeader));
00294 
00295     // Succeed, first record.
00296     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00297     validateRead1(samRecord);
00298 
00299     // Succeed, higher coordinate.
00300     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00301     validateRead2(samRecord);
00302 
00303     // Failed sort order - due to coord.
00304     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00305     validateRead3(samRecord);
00306 
00307     // Failed sort order - due to coord.
00308     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00309     validateRead4(samRecord);
00310 
00311     // Succeed, new reference id
00312     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00313     validateRead5(samRecord);
00314 
00315     // Fail, previous reference id.
00316     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00317     validateRead6(samRecord);
00318 
00319     // Succeed, same reference id, higher coord.
00320     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00321     validateRead7(samRecord);
00322 
00323     // Succeed, *, new reference id.
00324     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00325     validateRead8(samRecord);
00326 
00327     // Fail, reference id is not *
00328     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00329     validateRead9(samRecord);
00330 
00331     // Succeed, valid reference id, and no coordinate.
00332     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00333     validateRead10(samRecord);
00334 
00335 
00336     ////////////////////////////////////////////
00337     // Reopen the file for reading
00338     assert(inSam.OpenForRead("testFiles/testSamSOquery.sam"));
00339 
00340     // Set the validation to FLAG, SO set to queryname.
00341     inSam.setSortedValidation(SamFile::FLAG);
00342 
00343     // Read the SAM Header.
00344     assert(inSam.ReadHeader(samHeader));
00345    
00346     // Succeed, first record.
00347     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00348     validateRead1(samRecord);
00349 
00350     // Succeed, same name.
00351     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00352     validateRead2(samRecord);
00353 
00354     // Failed sort order - due to alphabetical
00355     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00356     validateRead3(samRecord);
00357 
00358     // Failed sort order - due to alphabetical
00359     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00360     validateRead4(samRecord);
00361 
00362     // Failed sort order - due to alphabetical
00363     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00364     validateRead5(samRecord);
00365 
00366     // Failed sort order - due to alphabetical
00367     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00368     validateRead6(samRecord);
00369 
00370     // Failed sort order - due to alphabetical
00371     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00372     validateRead7(samRecord);
00373 
00374     // Succeed
00375     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00376     validateRead8(samRecord);
00377 
00378     // Failed sort order - due to alphabetical
00379     assert(inSam.ReadRecord(samHeader, samRecord) == false);
00380     validateRead9(samRecord);
00381 
00382     // Succeed
00383     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00384     validateRead10(samRecord);
00385 
00386     ////////////////////////////////////////////
00387     // Reopen the file for reading, SO flag set to junk.
00388     assert(inSam.OpenForRead("testFiles/testSamSOinvalid.sam"));
00389 
00390     // Set the validation to the SO Flag.  Not set to anything valid,
00391     // so it is considered UNSORTED, so all reads should pass.
00392     inSam.setSortedValidation(SamFile::FLAG);
00393 
00394     // Read the SAM Header.
00395     assert(inSam.ReadHeader(samHeader));
00396    
00397     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00398     validateRead1(samRecord);
00399 
00400     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00401     validateRead2(samRecord);
00402 
00403     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00404     validateRead3(samRecord);
00405 
00406     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00407     validateRead4(samRecord);
00408 
00409     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00410     validateRead5(samRecord);
00411 
00412     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00413     validateRead6(samRecord);
00414 
00415     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00416     validateRead7(samRecord);
00417 
00418     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00419     validateRead8(samRecord);
00420 
00421     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00422     validateRead9(samRecord);
00423 
00424     assert(inSam.ReadRecord(samHeader, samRecord) == true);
00425     validateRead10(samRecord);
00426 }
00427 
00428 
00429 
00430 void validateRead1ModQuality(SamRecord& samRecord)
00431 {
00432     //////////////////////////////////////////
00433     // Validate Record 1
00434     // Create record structure for validating.
00435     int expectedBlockSize = 89;
00436     const char* expectedReferenceName = "1";
00437     const char* expectedMateReferenceName = "1";
00438     const char* expectedMateReferenceNameOrEqual = "=";
00439 
00440     bamRecordStruct* expectedRecordPtr =
00441         (bamRecordStruct *) malloc(expectedBlockSize + sizeof(int));
00442 
00443     char tag[3];
00444     char type;
00445     void* value;
00446     bamRecordStruct* bufferPtr;
00447     unsigned char* varPtr;
00448 
00449     expectedRecordPtr->myBlockSize = expectedBlockSize;
00450     expectedRecordPtr->myReferenceID = 0;
00451     expectedRecordPtr->myPosition = 1010;
00452     expectedRecordPtr->myReadNameLength = 23;
00453     expectedRecordPtr->myMapQuality = 0;
00454     expectedRecordPtr->myBin = 4681;
00455     expectedRecordPtr->myCigarLength = 2;
00456     expectedRecordPtr->myFlag = 73;
00457     expectedRecordPtr->myReadLength = 5;
00458     expectedRecordPtr->myMateReferenceID = 0;
00459     expectedRecordPtr->myMatePosition = 1010;
00460     expectedRecordPtr->myInsertSize = 0;
00461    
00462     // Check the alignment end
00463     assert(samRecord.get0BasedAlignmentEnd() == 1016);
00464     assert(samRecord.get1BasedAlignmentEnd() == 1017);
00465     assert(samRecord.getAlignmentLength() == 7);
00466     assert(samRecord.get0BasedUnclippedStart() == 1010);
00467     assert(samRecord.get1BasedUnclippedStart() == 1011);
00468     assert(samRecord.get0BasedUnclippedEnd() == 1016);
00469     assert(samRecord.get1BasedUnclippedEnd() == 1017);
00470 
00471     // Check the accessors.
00472     assert(samRecord.getBlockSize() == expectedRecordPtr->myBlockSize);
00473     assert(samRecord.getReferenceID() == expectedRecordPtr->myReferenceID);
00474     assert(strcmp(samRecord.getReferenceName(), expectedReferenceName) == 0);
00475     assert(samRecord.get1BasedPosition() == expectedRecordPtr->myPosition + 1);
00476     assert(samRecord.get0BasedPosition() == expectedRecordPtr->myPosition);
00477     assert(samRecord.getReadNameLength() == 
00478            expectedRecordPtr->myReadNameLength);
00479     assert(samRecord.getMapQuality() == expectedRecordPtr->myMapQuality);
00480     assert(samRecord.getBin() == expectedRecordPtr->myBin);
00481     assert(samRecord.getCigarLength() == expectedRecordPtr->myCigarLength);
00482     assert(samRecord.getFlag() == expectedRecordPtr->myFlag);
00483     assert(samRecord.getReadLength() == expectedRecordPtr->myReadLength);
00484     assert(samRecord.getMateReferenceID() ==
00485            expectedRecordPtr->myMateReferenceID);
00486     assert(strcmp(samRecord.getMateReferenceName(), 
00487                   expectedMateReferenceName) == 0);
00488     assert(strcmp(samRecord.getMateReferenceNameOrEqual(), 
00489                   expectedMateReferenceNameOrEqual) == 0);
00490     assert(samRecord.get1BasedMatePosition() == 
00491            expectedRecordPtr->myMatePosition + 1);
00492     assert(samRecord.get0BasedMatePosition() ==
00493            expectedRecordPtr->myMatePosition);
00494     assert(samRecord.getInsertSize() == expectedRecordPtr->myInsertSize);
00495     assert(strcmp(samRecord.getReadName(), "1:1011:F:255+17M15D20M") == 0);
00496     assert(strcmp(samRecord.getCigar(), "5M2D") == 0);
00497     assert(strcmp(samRecord.getSequence(), "CCGAA") == 0);
00498     assert(strcmp(samRecord.getQuality(), "ABCDE") == 0);
00499     assert(samRecord.getNumOverlaps(1010, 1017) == 5);
00500     assert(samRecord.getNumOverlaps(1010, 1016) == 5);
00501     assert(samRecord.getNumOverlaps(1012, 1017) == 3);
00502     assert(samRecord.getNumOverlaps(1015, 1017) == 0);
00503     assert(samRecord.getNumOverlaps(1017, 1010) == 0);
00504     assert(samRecord.getNumOverlaps(1013, 1011) == 0);
00505     assert(samRecord.getNumOverlaps(-1, 1017) == 5);
00506 
00507     // Reset the tag iter, since the tags have already been read.
00508     samRecord.resetTagIter();
00509 
00510     // Check the tags.
00511     assert(samRecord.getNextSamTag(tag, type, &value) == true);
00512     assert(tag[0] == 'A');
00513     assert(tag[1] == 'M');
00514     assert(type == 'i');
00515     assert(*(char*)value == 0);
00516     assert(samRecord.getNextSamTag(tag, type, &value) == true);
00517     assert(tag[0] == 'M');
00518     assert(tag[1] == 'D');
00519     assert(type == 'Z');
00520     assert(*(String*)value == "37");
00521     assert(samRecord.getNextSamTag(tag, type, &value) == true);
00522     assert(tag[0] == 'N');
00523     assert(tag[1] == 'M');
00524     assert(type == 'i');
00525     assert(*(char*)value == 0);
00526     assert(samRecord.getNextSamTag(tag, type, &value) == true);
00527     assert(tag[0] == 'X');
00528     assert(tag[1] == 'T');
00529     assert(type == 'A');
00530     assert(*(char*)value == 'R');
00531     // No more tags, should return false.
00532     assert(samRecord.getNextSamTag(tag, type, &value) == false);
00533     assert(samRecord.getNextSamTag(tag, type, &value) == false);
00534 
00535     // Get the record ptr.   
00536     bufferPtr = (bamRecordStruct*)samRecord.getRecordBuffer();
00537     // Validate the buffers match.
00538     assert(bufferPtr->myBlockSize == expectedRecordPtr->myBlockSize);
00539     assert(bufferPtr->myReferenceID == expectedRecordPtr->myReferenceID);
00540     assert(bufferPtr->myPosition == expectedRecordPtr->myPosition);
00541     assert(bufferPtr->myReadNameLength == expectedRecordPtr->myReadNameLength);
00542     assert(bufferPtr->myMapQuality == expectedRecordPtr->myMapQuality);
00543     assert(bufferPtr->myBin == expectedRecordPtr->myBin);
00544     assert(bufferPtr->myCigarLength == expectedRecordPtr->myCigarLength);
00545     assert(bufferPtr->myFlag == expectedRecordPtr->myFlag);
00546     assert(bufferPtr->myReadLength == expectedRecordPtr->myReadLength);
00547     assert(bufferPtr->myMateReferenceID ==
00548            expectedRecordPtr->myMateReferenceID);
00549     assert(bufferPtr->myMatePosition == expectedRecordPtr->myMatePosition);
00550     assert(bufferPtr->myInsertSize == expectedRecordPtr->myInsertSize);
00551 
00552     // Validate the variable length fields in the buffer.
00553     // Set the pointer to the start of the variable fields.
00554     varPtr = (unsigned char*)(&(bufferPtr->myData[0]));
00555 
00556     // Validate the readname.
00557     for(int i = 0; i < expectedRecordPtr->myReadNameLength; i++)
00558     {
00559         assert(*varPtr == samRecord.getReadName()[i]);
00560         varPtr++;
00561     }
00562 
00563     // Validate the cigar.
00564     // The First cigar is 5M which is 5 << 4 | 0 = 80
00565     assert(*(unsigned int*)varPtr == 80);
00566     // Increment the varptr the size of an int.
00567     varPtr += 4;
00568     // The 2nd cigar is 2D which is 2 << 4 | 2 = 34
00569     assert(*(unsigned int*)varPtr == 34);
00570     // Increment the varptr the size of an int.
00571     varPtr += 4;
00572    
00573     // Validate the sequence.
00574     // CC = 0x22
00575     assert(*varPtr == 0x22);
00576     varPtr++;
00577     // GA = 0x41
00578     assert(*varPtr == 0x41);
00579     varPtr++;
00580     // A  = 0x10
00581     assert(*varPtr == 0x10);
00582     varPtr++;
00583   
00584     // Validate the Quality
00585     for(int i = 0; i < expectedRecordPtr->myReadLength; i++)
00586     {
00587         assert(*varPtr == samRecord.getQuality()[i] - 33);
00588         varPtr++;
00589     }
00590 
00591     // Validate the tags.  
00592     assert(*varPtr == 'A');
00593     varPtr++;
00594     assert(*varPtr == 'M');
00595     varPtr++;
00596     assert(*varPtr == 'C');
00597     varPtr++;
00598     assert(*varPtr == 0);
00599     varPtr++;
00600     assert(*varPtr == 'M');
00601     varPtr++;
00602     assert(*varPtr == 'D');
00603     varPtr++;
00604     assert(*varPtr == 'Z');
00605     varPtr++;
00606     assert(*varPtr == '3');
00607     varPtr++;
00608     assert(*varPtr == '7');
00609     varPtr++;
00610     assert(*varPtr == 0);
00611     varPtr++;
00612     assert(*varPtr == 'N');
00613     varPtr++;
00614     assert(*varPtr == 'M');
00615     varPtr++;
00616     assert(*varPtr == 'C');
00617     varPtr++;
00618     assert(*varPtr == 0);
00619     varPtr++;
00620     assert(*varPtr == 'X');
00621     varPtr++;
00622     assert(*varPtr == 'T');
00623     varPtr++;
00624     assert(*varPtr == 'A');
00625     varPtr++;
00626     assert(*varPtr == 'R');
00627     varPtr++;
00628 }
00629 
00630 
00631 void testModHeader(SamFileHeader& samHeader)
00632 {
00633     // Check the header line.
00634     std::string headerString = "";
00635     assert(samHeader.getHeaderString(headerString) == true);
00636     assert(headerString == "@SQ\tSN:1\tLN:247249719\n@SQ\tSN:2\tLN:242951149\n@SQ\tSN:3\tLN:199501827\n@SQ\tSN:4\tLN:191273063\n@SQ\tSN:5\tLN:180857866\n@SQ\tSN:6\tLN:170899992\n@SQ\tSN:7\tLN:158821424\n@SQ\tSN:8\tLN:146274826\n@SQ\tSN:9\tLN:140273252\n@SQ\tSN:10\tLN:135374737\n@SQ\tSN:11\tLN:134452384\n@SQ\tSN:12\tLN:132349534\n@SQ\tSN:13\tLN:114142980\n@SQ\tSN:14\tLN:106368585\n@SQ\tSN:15\tLN:100338915\n@SQ\tSN:16\tLN:88827254\n@SQ\tSN:17\tLN:78774742\n@SQ\tSN:18\tLN:76117153\n@SQ\tSN:19\tLN:63811651\n@SQ\tSN:20\tLN:62435964\n@SQ\tSN:21\tLN:46944323\n@SQ\tSN:22\tLN:49691432\n@SQ\tSN:X\tLN:154913754\n@RG\tID:myID\tLB:library\tSM:sample\n@RG\tID:myID2\tSM:sample2\tLB:library2\n@CO\tComment 1\n@CO\tComment 2\n");
00637 
00638     // Remove a tag - by setting it to "".
00639     assert(samHeader.setRGTag("LB", "", "myID2") == true);
00640 
00641 
00642     // Check the header line.
00643     assert(samHeader.getHeaderString(headerString) == true);
00644     assert(headerString == "@SQ\tSN:1\tLN:247249719\n@SQ\tSN:2\tLN:242951149\n@SQ\tSN:3\tLN:199501827\n@SQ\tSN:4\tLN:191273063\n@SQ\tSN:5\tLN:180857866\n@SQ\tSN:6\tLN:170899992\n@SQ\tSN:7\tLN:158821424\n@SQ\tSN:8\tLN:146274826\n@SQ\tSN:9\tLN:140273252\n@SQ\tSN:10\tLN:135374737\n@SQ\tSN:11\tLN:134452384\n@SQ\tSN:12\tLN:132349534\n@SQ\tSN:13\tLN:114142980\n@SQ\tSN:14\tLN:106368585\n@SQ\tSN:15\tLN:100338915\n@SQ\tSN:16\tLN:88827254\n@SQ\tSN:17\tLN:78774742\n@SQ\tSN:18\tLN:76117153\n@SQ\tSN:19\tLN:63811651\n@SQ\tSN:20\tLN:62435964\n@SQ\tSN:21\tLN:46944323\n@SQ\tSN:22\tLN:49691432\n@SQ\tSN:X\tLN:154913754\n@RG\tID:myID\tLB:library\tSM:sample\n@RG\tID:myID2\tSM:sample2\n@CO\tComment 1\n@CO\tComment 2\n");
00645 
00646     //  Add an HD tag.
00647     SamHeaderHD* hd = new SamHeaderHD();
00648     assert(hd->setTag("VN", "1.3") == true);
00649     assert(samHeader.addHD(hd) == true);
00650     assert(strcmp(samHeader.getHDTagValue("VN"), "1.3") == 0);
00651     assert(samHeader.getHeaderString(headerString) == true);
00652     assert(headerString == "@SQ\tSN:1\tLN:247249719\n@SQ\tSN:2\tLN:242951149\n@SQ\tSN:3\tLN:199501827\n@SQ\tSN:4\tLN:191273063\n@SQ\tSN:5\tLN:180857866\n@SQ\tSN:6\tLN:170899992\n@SQ\tSN:7\tLN:158821424\n@SQ\tSN:8\tLN:146274826\n@SQ\tSN:9\tLN:140273252\n@SQ\tSN:10\tLN:135374737\n@SQ\tSN:11\tLN:134452384\n@SQ\tSN:12\tLN:132349534\n@SQ\tSN:13\tLN:114142980\n@SQ\tSN:14\tLN:106368585\n@SQ\tSN:15\tLN:100338915\n@SQ\tSN:16\tLN:88827254\n@SQ\tSN:17\tLN:78774742\n@SQ\tSN:18\tLN:76117153\n@SQ\tSN:19\tLN:63811651\n@SQ\tSN:20\tLN:62435964\n@SQ\tSN:21\tLN:46944323\n@SQ\tSN:22\tLN:49691432\n@SQ\tSN:X\tLN:154913754\n@RG\tID:myID\tLB:library\tSM:sample\n@RG\tID:myID2\tSM:sample2\n@HD\tVN:1.3\n@CO\tComment 1\n@CO\tComment 2\n");
00653 
00654     // Try adding another HD tag.
00655     SamHeaderHD* hd2 = new SamHeaderHD();
00656     assert(hd2->setTag("VN", "1.4") == true);
00657     assert(samHeader.addHD(hd2) == false);
00658     assert(strcmp(samHeader.getHDTagValue("VN"), "1.4") != 0);
00659     assert(strcmp(samHeader.getHDTagValue("VN"), "1.3") == 0);
00660     assert(samHeader.getHeaderString(headerString) == true);
00661     assert(headerString == "@SQ\tSN:1\tLN:247249719\n@SQ\tSN:2\tLN:242951149\n@SQ\tSN:3\tLN:199501827\n@SQ\tSN:4\tLN:191273063\n@SQ\tSN:5\tLN:180857866\n@SQ\tSN:6\tLN:170899992\n@SQ\tSN:7\tLN:158821424\n@SQ\tSN:8\tLN:146274826\n@SQ\tSN:9\tLN:140273252\n@SQ\tSN:10\tLN:135374737\n@SQ\tSN:11\tLN:134452384\n@SQ\tSN:12\tLN:132349534\n@SQ\tSN:13\tLN:114142980\n@SQ\tSN:14\tLN:106368585\n@SQ\tSN:15\tLN:100338915\n@SQ\tSN:16\tLN:88827254\n@SQ\tSN:17\tLN:78774742\n@SQ\tSN:18\tLN:76117153\n@SQ\tSN:19\tLN:63811651\n@SQ\tSN:20\tLN:62435964\n@SQ\tSN:21\tLN:46944323\n@SQ\tSN:22\tLN:49691432\n@SQ\tSN:X\tLN:154913754\n@RG\tID:myID\tLB:library\tSM:sample\n@RG\tID:myID2\tSM:sample2\n@HD\tVN:1.3\n@CO\tComment 1\n@CO\tComment 2\n");
00662 
00663     // Remove the entire HD Tag.
00664     assert(samHeader.removeHD() == true);
00665     assert(strcmp(samHeader.getHDTagValue("VN"), "") == 0);
00666     assert(samHeader.getHeaderString(headerString) == true);
00667     assert(headerString == "@SQ\tSN:1\tLN:247249719\n@SQ\tSN:2\tLN:242951149\n@SQ\tSN:3\tLN:199501827\n@SQ\tSN:4\tLN:191273063\n@SQ\tSN:5\tLN:180857866\n@SQ\tSN:6\tLN:170899992\n@SQ\tSN:7\tLN:158821424\n@SQ\tSN:8\tLN:146274826\n@SQ\tSN:9\tLN:140273252\n@SQ\tSN:10\tLN:135374737\n@SQ\tSN:11\tLN:134452384\n@SQ\tSN:12\tLN:132349534\n@SQ\tSN:13\tLN:114142980\n@SQ\tSN:14\tLN:106368585\n@SQ\tSN:15\tLN:100338915\n@SQ\tSN:16\tLN:88827254\n@SQ\tSN:17\tLN:78774742\n@SQ\tSN:18\tLN:76117153\n@SQ\tSN:19\tLN:63811651\n@SQ\tSN:20\tLN:62435964\n@SQ\tSN:21\tLN:46944323\n@SQ\tSN:22\tLN:49691432\n@SQ\tSN:X\tLN:154913754\n@RG\tID:myID\tLB:library\tSM:sample\n@RG\tID:myID2\tSM:sample2\n@CO\tComment 1\n@CO\tComment 2\n");
00668 
00669     // Remove an entire SQ Tag.
00670     assert(strcmp(samHeader.getSQTagValue("LN", "11"), "134452384") == 0);
00671     assert(samHeader.removeSQ("11") == true);
00672     assert(strcmp(samHeader.getSQTagValue("LN", "11"), "") == 0);
00673     assert(samHeader.getHeaderString(headerString) == true);
00674     assert(headerString == "@SQ\tSN:1\tLN:247249719\n@SQ\tSN:2\tLN:242951149\n@SQ\tSN:3\tLN:199501827\n@SQ\tSN:4\tLN:191273063\n@SQ\tSN:5\tLN:180857866\n@SQ\tSN:6\tLN:170899992\n@SQ\tSN:7\tLN:158821424\n@SQ\tSN:8\tLN:146274826\n@SQ\tSN:9\tLN:140273252\n@SQ\tSN:10\tLN:135374737\n@SQ\tSN:12\tLN:132349534\n@SQ\tSN:13\tLN:114142980\n@SQ\tSN:14\tLN:106368585\n@SQ\tSN:15\tLN:100338915\n@SQ\tSN:16\tLN:88827254\n@SQ\tSN:17\tLN:78774742\n@SQ\tSN:18\tLN:76117153\n@SQ\tSN:19\tLN:63811651\n@SQ\tSN:20\tLN:62435964\n@SQ\tSN:21\tLN:46944323\n@SQ\tSN:22\tLN:49691432\n@SQ\tSN:X\tLN:154913754\n@RG\tID:myID\tLB:library\tSM:sample\n@RG\tID:myID2\tSM:sample2\n@CO\tComment 1\n@CO\tComment 2\n");
00675 
00676     // Try adding a null HD tag.
00677     hd = NULL;
00678     assert(samHeader.addHD(hd) == false);
00679     assert(strcmp(samHeader.getHDTagValue("VN"), "") == 0);
00680     assert(strcmp(samHeader.getHDTagValue("VN"), "1.4") != 0);
00681     assert(strcmp(samHeader.getHDTagValue("VN"), "1.3") != 0);
00682     assert(samHeader.getHeaderString(headerString) == true);
00683     assert(headerString == "@SQ\tSN:1\tLN:247249719\n@SQ\tSN:2\tLN:242951149\n@SQ\tSN:3\tLN:199501827\n@SQ\tSN:4\tLN:191273063\n@SQ\tSN:5\tLN:180857866\n@SQ\tSN:6\tLN:170899992\n@SQ\tSN:7\tLN:158821424\n@SQ\tSN:8\tLN:146274826\n@SQ\tSN:9\tLN:140273252\n@SQ\tSN:10\tLN:135374737\n@SQ\tSN:12\tLN:132349534\n@SQ\tSN:13\tLN:114142980\n@SQ\tSN:14\tLN:106368585\n@SQ\tSN:15\tLN:100338915\n@SQ\tSN:16\tLN:88827254\n@SQ\tSN:17\tLN:78774742\n@SQ\tSN:18\tLN:76117153\n@SQ\tSN:19\tLN:63811651\n@SQ\tSN:20\tLN:62435964\n@SQ\tSN:21\tLN:46944323\n@SQ\tSN:22\tLN:49691432\n@SQ\tSN:X\tLN:154913754\n@RG\tID:myID\tLB:library\tSM:sample\n@RG\tID:myID2\tSM:sample2\n@CO\tComment 1\n@CO\tComment 2\n");
00684 
00685     // Try adding a null SQ tag.
00686     SamHeaderSQ* sq = NULL;
00687     assert(samHeader.addSQ(sq) == false);
00688     assert(samHeader.getHeaderString(headerString) == true);
00689     assert(headerString == "@SQ\tSN:1\tLN:247249719\n@SQ\tSN:2\tLN:242951149\n@SQ\tSN:3\tLN:199501827\n@SQ\tSN:4\tLN:191273063\n@SQ\tSN:5\tLN:180857866\n@SQ\tSN:6\tLN:170899992\n@SQ\tSN:7\tLN:158821424\n@SQ\tSN:8\tLN:146274826\n@SQ\tSN:9\tLN:140273252\n@SQ\tSN:10\tLN:135374737\n@SQ\tSN:12\tLN:132349534\n@SQ\tSN:13\tLN:114142980\n@SQ\tSN:14\tLN:106368585\n@SQ\tSN:15\tLN:100338915\n@SQ\tSN:16\tLN:88827254\n@SQ\tSN:17\tLN:78774742\n@SQ\tSN:18\tLN:76117153\n@SQ\tSN:19\tLN:63811651\n@SQ\tSN:20\tLN:62435964\n@SQ\tSN:21\tLN:46944323\n@SQ\tSN:22\tLN:49691432\n@SQ\tSN:X\tLN:154913754\n@RG\tID:myID\tLB:library\tSM:sample\n@RG\tID:myID2\tSM:sample2\n@CO\tComment 1\n@CO\tComment 2\n");
00690 
00691     // Try adding an HD tag again.
00692     assert(samHeader.addHD(hd2) == true);
00693     assert(strcmp(samHeader.getHDTagValue("VN"), "1.4") == 0);
00694     assert(strcmp(samHeader.getHDTagValue("VN"), "1.3") != 0);
00695     assert(samHeader.getHeaderString(headerString) == true);
00696     assert(headerString == "@SQ\tSN:1\tLN:247249719\n@SQ\tSN:2\tLN:242951149\n@SQ\tSN:3\tLN:199501827\n@SQ\tSN:4\tLN:191273063\n@SQ\tSN:5\tLN:180857866\n@SQ\tSN:6\tLN:170899992\n@SQ\tSN:7\tLN:158821424\n@SQ\tSN:8\tLN:146274826\n@SQ\tSN:9\tLN:140273252\n@SQ\tSN:10\tLN:135374737\n@SQ\tSN:12\tLN:132349534\n@SQ\tSN:13\tLN:114142980\n@SQ\tSN:14\tLN:106368585\n@SQ\tSN:15\tLN:100338915\n@SQ\tSN:16\tLN:88827254\n@SQ\tSN:17\tLN:78774742\n@SQ\tSN:18\tLN:76117153\n@SQ\tSN:19\tLN:63811651\n@SQ\tSN:20\tLN:62435964\n@SQ\tSN:21\tLN:46944323\n@SQ\tSN:22\tLN:49691432\n@SQ\tSN:X\tLN:154913754\n@RG\tID:myID\tLB:library\tSM:sample\n@RG\tID:myID2\tSM:sample2\n@HD\tVN:1.4\n@CO\tComment 1\n@CO\tComment 2\n");
00697 
00698 
00699     // TODO Get the comments.
00700 
00701 }
Generated on Tue Aug 23 18:19:04 2011 for libStatGen Software by  doxygen 1.6.3