Bam2Fastq - add option to get the quality from a tag

src/Bam2FastQ.cpp

@@ -40,6 +40,8 @@ Bam2FastQ::Bam2FastQ()
       myNumPairs(0),
       myNumUnpaired(0),
       mySplitRG(false),
+      myQField(""),
+      myNumQualTagErrors(0),
       myReverseComp(true),
       myRNPlus(false),
       myOutBase(""),
@@ -75,13 +77,15 @@ void Bam2FastQ::description()
 void Bam2FastQ::usage()
 {
     BamExecutable::usage();
-    std::cerr << "\t./bam bam2FastQ --in <inputFile> [--readName] [--refFile <referenceFile>] [--outBase <outputFileBase>] [--firstOut <1stReadInPairOutFile>] [--merge|--secondOut <2ndReadInPairOutFile>] [--unpairedOut <unpairedOutFile>] [--firstRNExt <firstInPairReadNameExt>] [--secondRNExt <secondInPairReadNameExt>] [--rnPlus] [--noReverseComp] [--noeof] [--params]" << std::endl;
+    std::cerr << "\t./bam bam2FastQ --in <inputFile> [--readName] [--splitRG] [--qualField <tag>] [--refFile <referenceFile>] [--outBase <outputFileBase>] [--firstOut <1stReadInPairOutFile>] [--merge|--secondOut <2ndReadInPairOutFile>] [--unpairedOut <unpairedOutFile>] [--firstRNExt <firstInPairReadNameExt>] [--secondRNExt <secondInPairReadNameExt>] [--rnPlus] [--noReverseComp] [--noeof] [--params]" << std::endl;
     std::cerr << "\tRequired Parameters:" << std::endl;
     std::cerr << "\t\t--in       : the SAM/BAM file to convert to FastQ" << std::endl;
     std::cerr << "\tOptional Parameters:" << std::endl;
     std::cerr << "\t\t--readname      : Process the BAM as readName sorted instead\n"
               << "\t\t                  of coordinate if the header does not indicate a sort order." << std::endl;
     std::cerr << "\t\t--splitRG       : Split into RG specific fastqs." << std::endl;
+    std::cerr << "\t\t--qualField     : Use the base quality from the specified tag\n";
+    std::cerr << "\t\t                  rather than from the Quality field (default)" << std::endl;
     std::cerr << "\t\t--merge         : Generate 1 interleaved (merged) FASTQ for paired-ends (unpaired in a separate file)\n"
               << "\t\t                  use firstOut to override the filename of the interleaved file." << std::endl;
     std::cerr << "\t\t--refFile       : Reference file for converting '=' in the sequence to the actual base" << std::endl;
@@ -127,6 +131,8 @@ int Bam2FastQ::execute(int argc, char **argv)
     myNumPairs = 0;
     myNumUnpaired = 0;
     mySplitRG = false;
+    myQField = "";
+    myNumQualTagErrors = 0;
     myReverseComp = true;
     myRNPlus = false;
     myFirstRNExt = DEFAULT_FIRST_EXT;
@@ -140,6 +146,7 @@ int Bam2FastQ::execute(int argc, char **argv)
         LONG_PARAMETER_GROUP("Optional Parameters")
         LONG_PARAMETER("readName", &readName)
         LONG_PARAMETER("splitRG", &mySplitRG)
+        LONG_STRINGPARAMETER("qualField", &myQField)
         LONG_PARAMETER("merge", &interleave)
         LONG_STRINGPARAMETER("refFile", &refFile)
         LONG_STRINGPARAMETER("firstRNExt", &myFirstRNExt)
@@ -403,6 +410,11 @@ int Bam2FastQ::execute(int argc, char **argv)
                   << " reads, so they were written as unpaired\n"
                   << "  (not included in either of the above counts).\n";
     }
+    if(myNumQualTagErrors != 0)
+    {
+        std::cerr << myNumQualTagErrors << " records did not have tag "
+                  << myQField.c_str() << " or it was invalid, so the quality field was used for those records.\n";
+    }
 
     return(returnStatus);
 }
@@ -623,7 +635,32 @@ void Bam2FastQ::writeFastQ(SamRecord& samRec, IFILE filePtr,
     flag = samRec.getFlag();
     const char* readName = samRec.getReadName();
     sequence = samRec.getSequence();
-    quality = samRec.getQuality();
+    if(myQField.IsEmpty())
+    {
+        // Read the quality from the quality field
+        quality = samRec.getQuality();
+    }
+    else
+    {
+        // Read Quality from the specified tag
+        const String* qTagPtr = samRec.getStringTag(myQField.c_str());
+        if((qTagPtr != NULL) && (qTagPtr->Length() == (int)sequence.length()))
+        {
+            // Use the tag value for quality
+            quality = qTagPtr->c_str();
+        }
+        else
+        {
+            // Tag was not found, so use the quality field.
+            ++myNumQualTagErrors;
+            if(myNumQualTagErrors == 1)
+            {
+                std::cerr << "Bam2FastQ: " << myQField.c_str() 
+                          << " tag was not found/invalid, so using the quality field in records without the tag\n";
+            }
+            quality = samRec.getQuality();
+        }
+    }
 
     if(SamFlag::isReverse(flag) && myReverseComp)
     {

src/Bam2FastQ.h

@@ -76,6 +76,8 @@ class Bam2FastQ : public BamExecutable
     int myNumUnpaired;
 
     bool mySplitRG;
+    String myQField;
+    uint64_t myNumQualTagErrors;
     bool myReverseComp;
     bool myRNPlus;
 

test/expected/testBam2FastQCoordFirstRGFail.log

@@ -1,13 +1,15 @@
 Version: 1.0.13; Built: Mon Jun  7 12:26:20 EDT 2015 by mktrost
 
  bam2FastQ - Convert the specified BAM file to fastQs.
-	./bam bam2FastQ --in <inputFile> [--readName] [--refFile <referenceFile>] [--outBase <outputFileBase>] [--firstOut <1stReadInPairOutFile>] [--merge|--secondOut <2ndReadInPairOutFile>] [--unpairedOut <unpairedOutFile>] [--firstRNExt <firstInPairReadNameExt>] [--secondRNExt <secondInPairReadNameExt>] [--rnPlus] [--noReverseComp] [--noeof] [--params]
+	./bam bam2FastQ --in <inputFile> [--readName] [--splitRG] [--qualField <tag>] [--refFile <referenceFile>] [--outBase <outputFileBase>] [--firstOut <1stReadInPairOutFile>] [--merge|--secondOut <2ndReadInPairOutFile>] [--unpairedOut <unpairedOutFile>] [--firstRNExt <firstInPairReadNameExt>] [--secondRNExt <secondInPairReadNameExt>] [--rnPlus] [--noReverseComp] [--noeof] [--params]
 	Required Parameters:
 		--in       : the SAM/BAM file to convert to FastQ
 	Optional Parameters:
 		--readname      : Process the BAM as readName sorted instead
 		                  of coordinate if the header does not indicate a sort order.
 		--splitRG       : Split into RG specific fastqs.
+		--qualField     : Use the base quality from the specified tag
+		                  rather than from the Quality field (default)
 		--merge         : Generate 1 interleaved (merged) FASTQ for paired-ends (unpaired in a separate file)
 		                  use firstOut to override the filename of the interleaved file.
 		--refFile       : Reference file for converting '=' in the sequence to the actual base
@@ -33,8 +35,8 @@ Version: 1.0.13; Built: Mon Jun  7 12:26:20 EDT 2015 by mktrost
 
 Input Parameters
          Required Parameters : --in [testFiles/testBam2FastQCoord.sam]
-         Optional Parameters : --readName, --splitRG [ON], --merge,
-                               --refFile [], --firstRNExt [/1],
+         Optional Parameters : --readName, --splitRG [ON], --qualField [],
+                               --merge, --refFile [], --firstRNExt [/1],
                                --secondRNExt [/2], --rnPlus,
                                --noReverseComp [ON], --gzip, --noeof, --params
    Optional OutputFile Names : --outBase [],

test/expected/testBam2FastQCoordRG.log

@@ -1,3 +1,4 @@
+Bam2FastQ: OQ tag was not found/invalid, so using the quality field in records without the tag
 Paired Read, 7 but couldn't find mate, so writing as unpaired (single-ended)
 Paired Read, 1 but couldn't find mate, so writing as unpaired (single-ended)
 Paired Read, 12 but couldn't find mate, so writing as unpaired (single-ended)
@@ -7,3 +8,4 @@ Found 20 read pairs.
 Found 4 unpaired reads.
 Failed to find mates for 3 reads, so they were written as unpaired
   (not included in either of the above counts).
+45 records did not have tag OQ or it was invalid, so the quality field was used for those records.

test/expected/testBam2FastQCoordRG.rg1_2.fastq

@@ -13,7 +13,7 @@ CCCTTTCCCAAAGGGAAATTTCCC
 @8/2
 CCCTTTCCCAAAGGGAAATTTCCC
 +
-,><><>><32908lkge<.<.523
+9><><>><32908lkge<.<.654
 @15/2
 CCCCTTTCCCAAAGGGAAATTCCC
 +
@@ -33,7 +33,7 @@ CCCCTTTCCCAAAGGGAAATTCCC
 @10/2
 CCCTTTCCCAAAGGGAAATTTCCC
 +
-,><><>><32908lkge<.<.523
+8><><>><32908lkge<.<.321
 @3/2
 CCCCTTTCCCAAAGGGAAATTCCC
 +

test/expected/testBam2FastQCoordRGgz.log

@@ -0,0 +1,9 @@
+Paired Read, 7 but couldn't find mate, so writing as unpaired (single-ended)
+Paired Read, 1 but couldn't find mate, so writing as unpaired (single-ended)
+Paired Read, 12 but couldn't find mate, so writing as unpaired (single-ended)
+Both reads of 26 are first fragment, so splitting one to be in the 2nd fastq.
+
+Found 20 read pairs.
+Found 4 unpaired reads.
+Failed to find mates for 3 reads, so they were written as unpaired
+  (not included in either of the above counts).

test/expected/testBam2FastQCoordRGgz.rg1_2.fastq

@@ -0,0 +1,56 @@
+@1/2
+CCCTTTCCCAAAGGGAAATTTCCC
++
+<.,><><>><32909lkgerw523
+@5/2
+CCCTTTCCCAAAGGGAAATTTCCC
++
+<.,><><>><32909lkgerw523
+@7/2
+CCCTTTCCCAAAGGGAAATTTCCC
++
+,><><>><32908lkge<.<.523
+@8/2
+CCCTTTCCCAAAGGGAAATTTCCC
++
+,><><>><32908lkge<.<.523
+@15/2
+CCCCTTTCCCAAAGGGAAATTCCC
++
+<.,><><>><329020kgerw523
+@2/2
+CCCCTTTCCCAAAGGGAAATTCCC
++
+<.,><><>><32909lkgerw523
+@11/2
+CCCCTTTCCCAAAGGGAAATTCCC
++
+<.,><><>><32909lkgerw523
+@4/2
+CCCCTTTCCCAAAGGGAAATTCCC
++
+<.,><><>><32909lkgerw523
+@10/2
+CCCTTTCCCAAAGGGAAATTTCCC
++
+,><><>><32908lkge<.<.523
+@3/2
+CCCCTTTCCCAAAGGGAAATTCCC
++
+<.,><><>><32909lkgerw523
+@13/2
+CCCCTTTCCCAAAGGGAAATTCCC
++
+<.,><><>><32909lkgerw523
+@30/2
+CCCTTTCCCAAAGGGAAATTTCCC
++
+<-,><><>><32909lkge<.523
+@25/2
+GGGAATTTCCCTTTGGGAAAGGGG
++
+325wregkl90923<>><><>,.<
+@26/2
+GGGAATTTCCCTTTGGGAAAGGGG
++
+325wregkl90923<>><><>,.<

test/expected/testBam2FastQCoordSecondRGFail.log

@@ -1,13 +1,15 @@
 Version: 1.0.13; Built: Mon Jun  7 12:26:20 EDT 2015 by mktrost
 
  bam2FastQ - Convert the specified BAM file to fastQs.
-	./bam bam2FastQ --in <inputFile> [--readName] [--refFile <referenceFile>] [--outBase <outputFileBase>] [--firstOut <1stReadInPairOutFile>] [--merge|--secondOut <2ndReadInPairOutFile>] [--unpairedOut <unpairedOutFile>] [--firstRNExt <firstInPairReadNameExt>] [--secondRNExt <secondInPairReadNameExt>] [--rnPlus] [--noReverseComp] [--noeof] [--params]
+	./bam bam2FastQ --in <inputFile> [--readName] [--splitRG] [--qualField <tag>] [--refFile <referenceFile>] [--outBase <outputFileBase>] [--firstOut <1stReadInPairOutFile>] [--merge|--secondOut <2ndReadInPairOutFile>] [--unpairedOut <unpairedOutFile>] [--firstRNExt <firstInPairReadNameExt>] [--secondRNExt <secondInPairReadNameExt>] [--rnPlus] [--noReverseComp] [--noeof] [--params]
 	Required Parameters:
 		--in       : the SAM/BAM file to convert to FastQ
 	Optional Parameters:
 		--readname      : Process the BAM as readName sorted instead
 		                  of coordinate if the header does not indicate a sort order.
 		--splitRG       : Split into RG specific fastqs.
+		--qualField     : Use the base quality from the specified tag
+		                  rather than from the Quality field (default)
 		--merge         : Generate 1 interleaved (merged) FASTQ for paired-ends (unpaired in a separate file)
 		                  use firstOut to override the filename of the interleaved file.
 		--refFile       : Reference file for converting '=' in the sequence to the actual base
@@ -33,8 +35,8 @@ Version: 1.0.13; Built: Mon Jun  7 12:26:20 EDT 2015 by mktrost
 
 Input Parameters
          Required Parameters : --in [testFiles/testBam2FastQCoord.sam]
-         Optional Parameters : --readName, --splitRG [ON], --merge,
-                               --refFile [], --firstRNExt [/1],
+         Optional Parameters : --readName, --splitRG [ON], --qualField [],
+                               --merge, --refFile [], --firstRNExt [/1],
                                --secondRNExt [/2], --rnPlus,
                                --noReverseComp [ON], --gzip, --noeof, --params
    Optional OutputFile Names : --outBase [], --firstOut [],

test/expected/testBam2FastQCoordUnpairRGFail.log

@@ -1,13 +1,15 @@
 Version: 1.0.13; Built: Mon Jun  7 12:26:20 EDT 2015 by mktrost
 
  bam2FastQ - Convert the specified BAM file to fastQs.
-	./bam bam2FastQ --in <inputFile> [--readName] [--refFile <referenceFile>] [--outBase <outputFileBase>] [--firstOut <1stReadInPairOutFile>] [--merge|--secondOut <2ndReadInPairOutFile>] [--unpairedOut <unpairedOutFile>] [--firstRNExt <firstInPairReadNameExt>] [--secondRNExt <secondInPairReadNameExt>] [--rnPlus] [--noReverseComp] [--noeof] [--params]
+	./bam bam2FastQ --in <inputFile> [--readName] [--splitRG] [--qualField <tag>] [--refFile <referenceFile>] [--outBase <outputFileBase>] [--firstOut <1stReadInPairOutFile>] [--merge|--secondOut <2ndReadInPairOutFile>] [--unpairedOut <unpairedOutFile>] [--firstRNExt <firstInPairReadNameExt>] [--secondRNExt <secondInPairReadNameExt>] [--rnPlus] [--noReverseComp] [--noeof] [--params]
 	Required Parameters:
 		--in       : the SAM/BAM file to convert to FastQ
 	Optional Parameters:
 		--readname      : Process the BAM as readName sorted instead
 		                  of coordinate if the header does not indicate a sort order.
 		--splitRG       : Split into RG specific fastqs.
+		--qualField     : Use the base quality from the specified tag
+		                  rather than from the Quality field (default)
 		--merge         : Generate 1 interleaved (merged) FASTQ for paired-ends (unpaired in a separate file)
 		                  use firstOut to override the filename of the interleaved file.
 		--refFile       : Reference file for converting '=' in the sequence to the actual base
@@ -33,8 +35,8 @@ Version: 1.0.13; Built: Mon Jun  7 12:26:20 EDT 2015 by mktrost
 
 Input Parameters
          Required Parameters : --in [testFiles/testBam2FastQCoord.sam]
-         Optional Parameters : --readName, --splitRG [ON], --merge,
-                               --refFile [], --firstRNExt [/1],
+         Optional Parameters : --readName, --splitRG [ON], --qualField [],
+                               --merge, --refFile [], --firstRNExt [/1],
                                --secondRNExt [/2], --rnPlus,
                                --noReverseComp [ON], --gzip, --noeof, --params
    Optional OutputFile Names : --outBase [], --firstOut [], --secondOut [],

test/expected/testBam2FastQReadNameRG.log

@@ -1,7 +1,9 @@
 Paired Read, 1 but couldn't find mate, so writing as unpaired (single-ended)
+Bam2FastQ: OQ tag was not found/invalid, so using the quality field in records without the tag
 Paired Read, 7 but couldn't find mate, so writing as unpaired (single-ended)
 
 Found 14 read pairs.
 Found 0 unpaired reads.
 Failed to find mates for 2 reads, so they were written as unpaired
   (not included in either of the above counts).
+28 records did not have tag OQ or it was invalid, so the quality field was used for those records.

test/expected/testBam2FastQReadNameRG.rg1_1.fastq

@@ -1,7 +1,7 @@
 @1/1
 GGGAAATTTCCCTTTGGGAAAGGG
 +
-325wregkl90923<>><><>,.<
+123.<.<egkl80923<>><><>8
 @2/1
 GGGAAATTTCCCTTTGGGAAAGGG
 +

test/expected/testBam2FastQReadNameRG.rg1_2.fastq

@@ -1,7 +1,7 @@
 @1/2
 CCCTTTCCCAAAGGGAAATTTCCC
 +
-<.,><><>><32909lkgerw523
+9><><>><32908lkge<.<.654
 @2/2
 CCCCTTTCCCAAAGGGAAATTCCC
 +

test/testBam2FastQ.sh

@@ -350,7 +350,7 @@ let "status |= $?"
 
 ##########################################
 # Test valid splitRG command - RG
-../bin/bam bam2FastQ --in testFiles/testBam2FastQCoordRG.sam --outBase results/testBam2FastQCoordRG --noph --splitRG 2> results/testBam2FastQCoordRG.log
+../bin/bam bam2FastQ --in testFiles/testBam2FastQCoordRG.sam --outBase results/testBam2FastQCoordRG --noph --splitRG --qualField OQ 2> results/testBam2FastQCoordRG.log
 let "status |= $?"
 diff results/testBam2FastQCoordRG_1.fastq expected/testBam2FastQCoordRG_1.fastq
 let "status |= $?"
@@ -379,7 +379,7 @@ then
   let "status = 1"
 fi
 
-../bin/bam bam2FastQ --readName --in testFiles/testBam2FastQReadNameRG.sam --outBase results/testBam2FastQReadNameRG --noph --splitRG 2> results/testBam2FastQReadNameRG.log
+../bin/bam bam2FastQ --readName --in testFiles/testBam2FastQReadNameRG.sam --outBase results/testBam2FastQReadNameRG --noph --splitRG --qualField OQ 2> results/testBam2FastQReadNameRG.log
 let "status |= $?"
 diff results/testBam2FastQReadNameRG_1.fastq expected/testBam2FastQReadNameRG_1.fastq
 let "status |= $?"
@@ -443,13 +443,13 @@ diff <(gunzip -c results/testBam2FastQCoordRGgz.rg1.fastq) expected/testBam2Fast
 let "status |= $?"
 diff <(gunzip -c results/testBam2FastQCoordRGgz.rg1_1.fastq) expected/testBam2FastQCoordRG.rg1_1.fastq
 let "status |= $?"
-diff <(gunzip -c results/testBam2FastQCoordRGgz.rg1_2.fastq) expected/testBam2FastQCoordRG.rg1_2.fastq
+diff <(gunzip -c results/testBam2FastQCoordRGgz.rg1_2.fastq) expected/testBam2FastQCoordRGgz.rg1_2.fastq
 let "status |= $?"
 diff <(gunzip -c results/testBam2FastQCoordRGgz.rg2_1.fastq) expected/testBam2FastQCoordRG.rg2_1.fastq
 let "status |= $?"
 diff <(gunzip -c results/testBam2FastQCoordRGgz.rg2_2.fastq) expected/testBam2FastQCoordRG.rg2_2.fastq
 let "status |= $?"
-diff results/testBam2FastQCoordRGgz.log expected/testBam2FastQCoordRG.log
+diff results/testBam2FastQCoordRGgz.log expected/testBam2FastQCoordRGgz.log
 let "status |= $?"
 diff results/testBam2FastQCoordRGgz.list expected/testBam2FastQCoordRGgz.list
 let "status |= $?"

test/testFiles/testBam2FastQCoordRG.sam

@@ -4,8 +4,8 @@
 @SQ	SN:4	LN:191273063
 @RG	ID:rg1	SM:sm1
 @RG	ID:rg2	PL:ILLUMINA	LB:lib2	SM:sm2	CN:cn2
-10	147	1	8	0	3H3S3M3D3M3I3M3P3M3D3M3S3H	=	100	0	gggAAATTTCCCTTTGGGAAAggg	325.<.<egkl80923<>><><>,	XZ:Z:forward/reverseEntireClip	RG:Z:rg1
-8	147	1	8	0	3H3S3M3D3M3I3M3P3M3D3M3S3H	=	10	0	gggAAATTTCCCTTTGGGAAAggg	325.<.<egkl80923<>><><>,	XZ:Z:forward/reverseClip	RG:Z:rg1
+10	147	1	8	0	3H3S3M3D3M3I3M3P3M3D3M3S3H	=	100	0	gggAAATTTCCCTTTGGGAAAggg	325.<.<egkl80923<>><><>,	XZ:Z:forward/reverseEntireClip	RG:Z:rg1	OQ:Z:123.<.<egkl80923<>><><>8
+8	147	1	8	0	3H3S3M3D3M3I3M3P3M3D3M3S3H	=	10	0	gggAAATTTCCCTTTGGGAAAggg	325.<.<egkl80923<>><><>,	XZ:Z:forward/reverseClip	RG:Z:rg1	OQ:Z:456.<.<egkl80923<>><><>9
 9	147	1	8	0	3H3S3M3D3M3I3M3P3M3D3M3S3H	=	10	0	gggAAATTTCCCTTTGGGAAAggg	325.<.<egkl80923<>><><>,	XZ:Z:forward/reverseClip	RG:Z:rg2
 7	147	1	8	0	3H3S3M3D3M3I3M3P3M3D3M3S3H	=	10	0	gggAAATTTCCCTTTGGGAAAggg	325.<.<egkl80923<>><><>,	XZ:Z:forward/reverseClip	RG:Z:rg1
 1	99	1	10	0	3H3S3M3D3M3I3M3P3M3D3M3S3H	=	10	0	gggAAATTTCCCTTTGGGAAAggg	325wregkl90923<>><><>,.<	XZ:Z:CompleteOverlap	RG:Z:rg1

test/testFiles/testBam2FastQReadNameRG.sam

@@ -4,8 +4,8 @@
 @SQ	SN:4	LN:191273063
 @RG	ID:rg1	SM:sm1
 @RG	ID:rg2	PL:ILLUMINA	LB:lib2	SM:sm2	CN:cn2
-1	99	1	10	0	3H3S3M3D3M3I3M3P3M3D3M3S3H	=	10	0	gggAAATTTCCCTTTGGGAAAggg	325wregkl90923<>><><>,.<	XZ:Z:CompleteOverlap	RG:Z:rg1
-1	147	1	10	0	3H3S3M3D3M3I3M3P3M3D3M3S3H	=	10	0	gggAAATTTCCCTTTGGGAAAggg	325wregkl90923<>><><>,.<	XZ:Z:CompleteOverlap	RG:Z:rg1
+1	99	1	10	0	3H3S3M3D3M3I3M3P3M3D3M3S3H	=	10	0	gggAAATTTCCCTTTGGGAAAggg	325wregkl90923<>><><>,.<	XZ:Z:CompleteOverlap	RG:Z:rg1	OQ:Z:123.<.<egkl80923<>><><>8
+1	147	1	10	0	3H3S3M3D3M3I3M3P3M3D3M3S3H	=	10	0	gggAAATTTCCCTTTGGGAAAggg	325wregkl90923<>><><>,.<	XZ:Z:CompleteOverlap	RG:Z:rg1	OQ:Z:456.<.<egkl80923<>><><>9
 1	147	1	10	0	3H3S3M3D3M3I3M3P3M3D3M3S3H	=	10	0	gggAAATTTCCCTTTGGGAAAggg	325wregkl90923<>><><>,.<	XZ:Z:CompleteOverlap No Change, no match	RG:Z:rg1
 2	99	1	10	0	3H3S3M3D3M3I3M3P3M3D3M3S3H	=	11	0	gggAAATTTCCCTTTGGGAAAggg	325wregkl90923<>><><>,.<	XZ:Z:Shift1	RG:Z:rg1
 2	147	1	11	0	3H3S2M3D3M3I3M3P3M3D4M3S3H	=	10	0	gggAATTTCCCTTTGGGAAAGggg	325wregkl90923<>><><>,.<	XZ:Z:Shift1	RG:Z:rg1