module Bio.Sequence.FastQ
(
readFastQ, hReadFastQ, parse
, writeFastQ, hWriteFastQ, unparse
, readSangerQ, hReadSangerQ
, writeSangerQ, hWriteSangerQ
, readIllumina, hReadIllumina
, writeIllumina, hWriteIllumina
) where
import System.IO
import qualified Data.ByteString.Lazy.Char8 as B
import qualified Data.ByteString.Lazy as BB
import Data.List (unfoldr)
import Bio.Sequence.SeqData
readSangerQ, readIllumina :: FilePath -> IO [Sequence Nuc]
readSangerQ = readFastQ
readIllumina f = addQual (negate 31) `fmap` readFastQ f
hReadSangerQ, hReadIllumina :: Handle -> IO [Sequence Nuc]
hReadSangerQ = hReadFastQ
hReadIllumina h = addQual (negate 31) `fmap` hReadFastQ h
writeSangerQ, writeIllumina :: FilePath -> [Sequence Nuc] -> IO ()
writeSangerQ = writeFastQ
writeIllumina f = writeFastQ f . addQual 31
hWriteSangerQ, hWriteIllumina :: Handle -> [Sequence Nuc] -> IO ()
hWriteSangerQ = hWriteFastQ
hWriteIllumina h = hWriteFastQ h . addQual 31
addQual :: Qual -> [Sequence Nuc] -> [Sequence Nuc]
addQual q = map (\(Seq h d mq) -> (Seq h d $ BB.map (+q) `fmap` mq))
readFastQ :: FilePath -> IO [Sequence Nuc]
readFastQ f = (go . B.lines) `fmap` B.readFile f
hReadFastQ :: Handle -> IO [Sequence Nuc]
hReadFastQ h = (go . B.lines) `fmap` B.hGetContents h
go :: [B.ByteString] -> [Sequence Nuc]
go = map (either error id) . unfoldr parse
parse :: [B.ByteString] -> Maybe (Either String (Sequence Nuc), [B.ByteString])
parse (h1:sd:h2:sq:rest) =
case (B.uncons h1,B.uncons h2) of
(Just ('@',h1name), Just ('+',h2name))
| h1name == h2name || B.null h2name
-> Just (Right $ Seq h1name sd (Just (BB.map (subtract 33) sq)), rest)
| otherwise
-> Just (Left $ "Bio.Sequence.FastQ: name mismatch:" ++ showStanza, rest)
_ -> Just (Left $ "Bio.Sequence.FastQ: illegal FastQ format:" ++ showStanza, rest)
where showStanza = unlines $ map B.unpack [ h1, sd, h2, sq ]
parse [] = Nothing
parse fs = let showStanza = unlines (map B.unpack fs)
err = Left $ "Bio.Sequence.FastQ: illegal number of lines in FastQ format: " ++ showStanza
in Just (err, [])
writeFastQ :: FilePath -> [Sequence Nuc] -> IO ()
writeFastQ f = B.writeFile f . B.concat . map unparse
hWriteFastQ :: Handle -> [Sequence Nuc] -> IO ()
hWriteFastQ h = B.hPut h . B.concat . map unparse
unparse :: Sequence Nuc -> B.ByteString
unparse (Seq h sd (Just sq)) =
B.unlines [B.cons '@' h, sd, B.cons '+' h, BB.map (+33) sq]
unparse (Seq h _ Nothing) = error ("Bio.Sequence.FastQ: sequence " ++ show (B.unpack h)
++" doesn not have quality data!")