libStatGen Software
1
|
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 "SamFile.h" 00019 #include "SamFlag.h" 00020 #include "Modify.h" 00021 00022 void testModify() 00023 { 00024 modify modTest; 00025 modTest.testModify("testFiles/testSam.sam"); 00026 #ifdef __ZLIB_AVAILABLE__ 00027 modTest.testModify("testFiles/testBam.bam"); 00028 #endif 00029 } 00030 00031 00032 void modify::testModify(const char* filename) 00033 { 00034 myFilename = filename; 00035 00036 modifyPosition(); 00037 modifyCigar(); 00038 00039 modifyFlag(); 00040 00041 modifyTags(); 00042 } 00043 00044 void modify::modifyPosition() 00045 { 00046 openAndRead1Rec(); 00047 00048 // Verify the initial bin. 00049 assert(samRecord.getBin() == 4681); 00050 00051 // Change the position and verify that the bin is updated. 00052 assert(samRecord.set0BasedPosition(33768)); 00053 00054 // Verify the bin was updated. 00055 assert(samRecord.getBin() == 4683); 00056 assert(samRecord.get0BasedPosition() == 33768); 00057 } 00058 00059 00060 void modify::modifyCigar() 00061 { 00062 openAndRead1Rec(); 00063 00064 // Verify the initial bin. 00065 assert(samRecord.getBin() == 4681); 00066 00067 // Change the Cigar such that it modifies the bin. 00068 assert(samRecord.setCigar("33768M")); 00069 00070 // Verify the bin was updated. 00071 assert(samRecord.getBin() == 585); 00072 } 00073 00074 00075 void modify::modifyFlag() 00076 { 00077 openAndRead1Rec(); 00078 00079 // Verify the initial bin. 00080 uint16_t flag = 73; 00081 assert(samRecord.getFlag() == flag); 00082 00083 SamFlag::setDuplicate(flag); 00084 assert(flag == 1097); 00085 assert(samRecord.setFlag(flag)); 00086 assert(samRecord.getFlag() == 1097); 00087 00088 SamFlag::setNotDuplicate(flag); 00089 assert(flag == 73); 00090 assert(samRecord.setFlag(flag)); 00091 assert(samRecord.getFlag() == 73); 00092 } 00093 00094 00095 void modify::openAndRead1Rec() 00096 { 00097 // Open the file for reading. 00098 assert(samIn.OpenForRead(myFilename.c_str())); 00099 00100 // Read the sam header. 00101 assert(samIn.ReadHeader(samHeader)); 00102 00103 // Read the first record. 00104 assert(samIn.ReadRecord(samHeader, samRecord)); 00105 } 00106 00107 00108 void modify::modifyTags() 00109 { 00110 assert(samIn.OpenForRead(myFilename.c_str())); 00111 // Read the sam header. 00112 assert(samIn.ReadHeader(samHeader)); 00113 00114 SamFile samOut; 00115 SamFile bamOut; 00116 00117 std::string inputType = myFilename.substr(myFilename.find_last_of('.')); 00118 std::string outFileBase = "results/updateTagFrom"; 00119 if(inputType == ".bam") 00120 { 00121 outFileBase += "Bam"; 00122 } 00123 else 00124 { 00125 outFileBase += "Sam"; 00126 } 00127 00128 std::string outFile = outFileBase + ".sam"; 00129 assert(samOut.OpenForWrite(outFile.c_str())); 00130 outFile = outFileBase + ".bam"; 00131 assert(bamOut.OpenForWrite(outFile.c_str())); 00132 assert(samOut.WriteHeader(samHeader)); 00133 assert(bamOut.WriteHeader(samHeader)); 00134 00135 int count = 0; 00136 // Read the records. 00137 while(samIn.ReadRecord(samHeader, samRecord)) 00138 { 00139 if(count == 0) 00140 { 00141 assert(samRecord.rmTag("MD", 'Z')); 00142 } 00143 else if(count == 2) 00144 { 00145 assert(samRecord.rmTags("XT:A;MD:Z;AB:c;NM:i")); 00146 } 00147 else if(count == 4) 00148 { 00149 assert(samRecord.rmTags("MD:Z,AB:c,NM:i")); 00150 } 00151 00152 assert(bamOut.WriteRecord(samHeader, samRecord)); 00153 assert(samOut.WriteRecord(samHeader, samRecord)); 00154 ++count; 00155 } 00156 }