diff --git a/biophd.cabal b/biophd.cabal --- a/biophd.cabal +++ b/biophd.cabal @@ -1,23 +1,23 @@ Name: biophd -Version: 0.0.5 +Version: 0.0.6 Synopsis: Library for reading phd sequence files Description: Library for reading phd sequence files -Homepage: https://patch-tag.com/r/dfornika/biophd/home +Homepage: https://github.com/dfornika/biophd/wiki License: GPL License-file: LICENSE Cabal-Version: >=1.6 -Author: Dan Fornika -Maintainer: dfornika@gmail.com +Author: Ketil Malde +Maintainer: Dan Fornika Stability: Provisional Category: Bioinformatics Build-Type: Simple Library - Build-depends: base >= 2 && < 5, biocore, bytestring, parsec, text, binary - Exposed-modules: Bio.Sequence.Phd, Bio.Sequence.PhdData, Bio.Sequence.PhdTag + Build-depends: base >= 2 && < 5, biocore, bytestring, parsec, text, binary, time, old-locale + Exposed-modules: Bio.Sequence.Phd Hs-source-dirs: src source-Repository head - type: darcs - location: http://www.patch-tag.com/r/dfornika/biophd + type: git + location: https://github.com/dfornika/biophd.git diff --git a/src/Bio/Sequence/Phd.hs b/src/Bio/Sequence/Phd.hs --- a/src/Bio/Sequence/Phd.hs +++ b/src/Bio/Sequence/Phd.hs @@ -1,29 +1,36 @@ -module Bio.Sequence.Phd(Phd(..), readPhd, readPhdTags) where +module Bio.Sequence.Phd( + Phd(..), + readPhd, hReadPhd, + readPhdTags) where import Bio.Core.Sequence import Bio.Sequence.PhdData -import qualified Bio.Sequence.PhdTag as PT +import Bio.Sequence.PhdTag import Text.ParserCombinators.Parsec hiding (label) -import qualified Data.ByteString as B -import qualified Data.ByteString.Lazy as LB -import qualified Data.ByteString.Lazy.Char8 as LBC +import qualified Data.ByteString as BBB +import qualified Data.ByteString.Lazy as BB +import qualified Data.ByteString.Lazy.Char8 as B import Data.Ix import Data.Int (Int64) +import Data.Char (isSpace) import Data.List import Data.Maybe import Data.Text (Text) +import Data.Time +import Data.Time.Format import Data.Binary (encode) import System.IO +import System.Locale -- | Parse a .phd file, extracting the contents as a PHD readPhd :: FilePath -> IO Phd readPhd f = return . mkPhd =<< readFile f -readPhdTags :: FilePath -> IO (Maybe [PT.PhdTag]) -readPhdTags f = return . (mkPhdTags . lines) =<< readFile f +readPhdTags :: FilePath -> IO (Maybe [PhdTag]) +readPhdTags f = return . mkPhdTags . lines =<< readFile f -- | Parse .phd contents from a handle hReadPhd :: Handle -> IO Phd @@ -32,51 +39,57 @@ -- | The actual phd parser. mkPhd :: String -> Phd -mkPhd inp = +mkPhd inp = let (hd:fs) = filter (not . null) . lines $ inp (magic,label) = splitAt 15 hd (comment'',seqAndTags) = break (== "BEGIN_DNA") fs comment' = init $ tail comment'' - comment = map (tail . snd) (map (break (==' ')) comment') + comment = map ((tail . snd) . break (== ' ')) comment' (seq', tags') = break (== "END_DNA") seqAndTags seq = drop 1 seq' tags = init $ tail tags' fields = words . unlines $ comment sdata = filter ((==3).length) . map words $ seq err = error "failed to parse quality value" - dna = concat $ map (!!0) $ sdata - qual = map ((\x -> read x :: Int) . (!!1)) $ sdata - traceInd = map ((\x -> read x :: Int) . (!!2)) $ sdata - in if (magic == "BEGIN_SEQUENCE ") then Phd (mkComment comment) - (mkDNABlock label dna qual traceInd) - (mkPhdTags tags) - - else error "Incorrectly formatted PHD file - missing BEGIN_SEQUENCE" + dna = concatMap (!! 0) sdata + qual = map ((\x -> read x :: Int) . (!! 1)) sdata + traceInd = map ((\x -> read x :: Int) . (!! 2)) sdata + in if magic == "BEGIN_SEQUENCE " then + Phd (mkComment comment) + (mkDNABlock label dna qual traceInd) + (mkPhdTags tags) + else + error "Incorrectly formatted PHD file - missing BEGIN_SEQUENCE" -- Todo: also check that we have a BEGIN_DNA/END_DNA region there. mkComment :: [String] -> Comment -mkComment com = Comment { chromatFile = com!!0 +mkComment com = Comment { chromatFile = head com , abiThumbprint = com!!1 , phredVersion = com!!2 , callMethod = com!!3 , qualityLevels = read (com!!4) :: Int - , time = com!!5 + , time = readTime defaultTimeLocale "%a %b %-e %T %Y" $ com!!5 , traceArrayMinIndex = read (com!!6) :: Int , traceArrayMaxIndex = read (com!!7) :: Int - , trim = if (length com > 10) then Just (com!!8) else Nothing - , chem = if (length com > 10) then (com!!9) else (com!!8) - , dye = if (length com > 10) then (com!!10) else (com !!9) } + , trim = if length com > 10 then + Just $ com!!8 + else + Nothing + , chem = (!!) com (if length com > 10 then 9 else 8) + , dye = (!!) com (if length com > 10 then 10 else 9) + } mkDNABlock :: String -> String -> [Int] -> [Int] -> DNABlock -mkDNABlock l d q t = DNABlock { label = l - , bases = SeqData { unSD = LBC.pack d } - , qualities = QualData { unQD = encode q } +mkDNABlock l d q t = DNABlock { label = SeqLabel { unSL = B.pack l} + , bases = SeqData { unSD = B.pack d } + , qualities = QualData { unQD = encode q } , traceIndices = t } -mkPhdTags :: [String] -> Maybe [PT.PhdTag] +mkPhdTags :: [String] -> Maybe [PhdTag] mkPhdTags phdLines = case groupByTags phdLines of [] -> Nothing - _ -> Just (map (fromJust . mkOnePhdTag) (groupByTags phdLines)) + _ -> Just (map (fromJust . mkOnePhdTag) + (groupByTags phdLines)) groupByTags :: [String] -> [[String]] groupByTags [] = [] @@ -84,22 +97,29 @@ let begin_indices = elemIndices "BEGIN_TAG" xs end_indices = elemIndices "END_TAG" xs tag_spans = zip begin_indices end_indices - grouping = \x -> drop (fst (tag_spans!!x)) (take (snd (tag_spans!!x) + 1) xs) - in map grouping (Data.Ix.range (0, (length tag_spans) -1)) + grouping x = drop (fst (tag_spans !! x)) (take (snd (tag_spans !! x) + 1) xs) + in map grouping (Data.Ix.range (0, length tag_spans - 1)) -mkOnePhdTag :: [String] -> Maybe PT.PhdTag -mkOnePhdTag td = case length td of +parseReadPosition :: String -> (Offset, Offset) +parseReadPosition line = (head positions, (head . tail) positions) + where positions = map (\x -> Offset {unOff = read x :: Int64}) $ dropLabel line + +dropLabel :: String -> [String] +dropLabel = tail . words + +mkOnePhdTag :: [String] -> Maybe PhdTag +mkOnePhdTag td = case length td of 0 -> Nothing - 9 -> Just PT.PhdTag { PT.tagType = drop 6 (td!!1) - , PT.source = drop 8 (td!!2) - , PT.unpaddedReadPosition = map (\x -> Offset {unOff = read x :: Int64}) (words (drop 19 (td!!3))) - , PT.date = drop 6 (td!!4) - , PT.comment = if td!!6 == "BEGIN_TAG" then "" - else td!!6 - } - _ -> Just PT.PhdTag { PT.tagType = drop 6 (td!!1) - , PT.source = drop 8 (td!!2) - , PT.unpaddedReadPosition = map (\x -> Offset {unOff = read x :: Int64}) (words (drop 19 (td!!3))) - , PT.date = drop 6 (td!!4) - , PT.comment = "" - } + 9 -> Just PhdTag { tagType = head $ dropLabel (td!!1) + , source = head $ dropLabel (td!!2) + , unpaddedReadPosition = parseReadPosition (td!!3) + , date = readTime defaultTimeLocale "%y/%d/%m %T" $ intercalate " " $ dropLabel (td!!4) + , Bio.Sequence.PhdTag.comment = if td!!6 == "BEGIN_TAG" then Nothing + else Just (td!!6) + } + _ -> Just PhdTag { tagType = head $ dropLabel (td!!1) + , source = head $ dropLabel (td!!2) + , unpaddedReadPosition = parseReadPosition (td!!3) + , date = readTime defaultTimeLocale "%y/%d/%m %T" $ intercalate " " $ dropLabel (td!!4) + , Bio.Sequence.PhdTag.comment = Nothing + } diff --git a/src/Bio/Sequence/PhdData.hs b/src/Bio/Sequence/PhdData.hs deleted file mode 100644 --- a/src/Bio/Sequence/PhdData.hs +++ /dev/null @@ -1,89 +0,0 @@ -module Bio.Sequence.PhdData where - -import Data.Maybe(fromJust) -import Bio.Core.Sequence -import qualified Bio.Sequence.PhdTag as PT -import qualified Data.ByteString.Lazy as LB -import qualified Data.ByteString.Lazy.Char8 as LBC - -{-- A .phd file consists of a DNA block with base and quality - values, followed by one or more (optional) tag blocks. --} -data Phd = Phd { comment :: Comment - , dnaBlock :: DNABlock - , phdTags :: Maybe [PT.PhdTag] - } deriving (Show) - -{-- These types are subject to change if it improves functionality, - but for now it's simplest to just call them String, Int etc.--} -data Comment = Comment - { chromatFile :: FilePath - , abiThumbprint :: String - , phredVersion :: String - , callMethod :: String - , qualityLevels :: Int - , time :: String - , traceArrayMinIndex :: Int - , traceArrayMaxIndex :: Int - , trim :: Maybe String - , chem :: String - , dye :: String - } deriving (Eq) - -instance Show Comment where - show (Comment cf abit pv cm ql ti mintai maxtai tr ch dy) = - ("\n" ++) $ unlines - [ "CHROMAT_FILE: " ++ cf - , "ABI_THUMBPRINT: " ++ abit - , "PHRED_VERSION: " ++ pv - , "CALL_METHOD: " ++ cm - , "QUALITY_LEVELS: " ++ show ql - , "TIME: " ++ show ti - , "TRACE_ARRAY_MIN_INDEX: " ++ show mintai - , "TRACE_ARRAY_MAX_INDEX: " ++ show maxtai - , "TRIM: " ++ if (tr == Nothing) then "" else fromJust tr - , "CHEM: " ++ ch - , "DYE: " ++ dy - ] - -data DNABlock = DNABlock - { label :: String - , bases :: SeqData - , qualities :: QualData - , traceIndices :: [Int] - } - -instance Show DNABlock where - show = LBC.unpack . toFasta - -instance BioSeq DNABlock where - seqlabel db = SeqLabel $ LBC.pack $ label db - seqdata db = bases db - seqlength db = Offset $ LBC.length $ unSD $ bases db - -instance BioSeqQual DNABlock where - seqqual = qualities - --- Some default values for the data types, useful for debugging in ghci - -defaultComment = Comment { chromatFile = "" - , abiThumbprint = "0" - , phredVersion = "0.980904.e" - , callMethod = "phred" - , qualityLevels = 99 - , time = "" - , traceArrayMinIndex = 0 - , traceArrayMaxIndex = 1 - , trim = Nothing - , chem = "unknown" - , dye = "unknown" } - -defaultDNABlock = DNABlock { label = "some_dna" - , bases = SeqData $ LBC.pack "aatgcatcta" - , qualities = QualData $ LBC.pack "0000000000" - , traceIndices = [0,1,2,3,4,5,6,7,8,9,10] } - -defaultPhdTag = PT.PhdTag { PT.tagType = "polymorphism" - , PT.source = "polyphred" - , PT.unpaddedReadPosition = [5, 5] - , PT.date = "01/01/70 00:00:00" - , PT.comment = "" } diff --git a/src/Bio/Sequence/PhdTag.hs b/src/Bio/Sequence/PhdTag.hs deleted file mode 100644 --- a/src/Bio/Sequence/PhdTag.hs +++ /dev/null @@ -1,19 +0,0 @@ -module Bio.Sequence.PhdTag where - -import Bio.Core.Sequence (Offset, unOff) - -data PhdTag = PhdTag - { tagType :: String - , source :: String - , unpaddedReadPosition :: [Offset] - , date :: String - , comment :: String - } deriving (Eq) - -instance Show PhdTag where - show (PhdTag tt so urp da co) = ("\n" ++) $ unlines $ map (" " ++) - [ "TYPE: " ++ show tt - , "SOURCE: " ++ show so - , "UNPADDED_READ_POS: " ++ show (map unOff urp) - , "DATE: " ++ show da - , "COMMENT: " ++ show co ]