BiobaseInfernal-0.8.1.0: Biobase/SElab/Model/Import.hs
-- | Import SElab HMM and CM models, and combine these together into full
-- CM models.
--
-- TODO A pair @Right CM +++ Left HMM@ should become a @Right CM@, a lonely
-- @Left HMM@ should become a @Left HMM@, a lonely @Right CM@ should become
-- a @Right CM@.
module Biobase.SElab.Model.Import where
import Control.Applicative
import Control.DeepSeq
import Control.Lens (set,over,(^.),zoom)
import Control.Monad (liftM2)
import Control.Monad.Trans.Class (lift)
import Control.Monad.Trans.Writer.Strict
import Control.Monad (void,unless)
import Control.Monad (when,replicateM)
import Control.Parallel.Strategies (using,parList,rdeepseq,parMap)
import Data.ByteString (ByteString)
import Data.Maybe (catMaybes)
import Data.Monoid
import Data.String (IsString)
import Data.Text (Text)
import Data.Text.Encoding (decodeUtf8)
import Debug.Trace
import qualified Data.Attoparsec.ByteString.Char8 as ABC
import qualified Data.Attoparsec.Text as AT
import qualified Data.ByteString.Char8 as BS
import qualified Data.List as L
import qualified Data.Map.Strict as M
import qualified Data.Text as T
import qualified Data.Text.IO as T
import qualified Pipes as P
import qualified Pipes.Attoparsec as PA
import qualified Pipes.ByteString as PB
import qualified Pipes.GZip as PG
import qualified Pipes.Parse as PP
import qualified Pipes.Prelude as P
import qualified Pipes.Safe.Prelude as PSP
import System.FilePath (takeExtension)
import System.IO (stdin,withFile,IOMode(..))
import Biobase.Types.Accession
import Pipes.Split.ByteString
import Biobase.SElab.CM.Import as CM
import Biobase.SElab.CM.Types as CM
import Biobase.SElab.HMM.Import as HMM
import Biobase.SElab.HMM.Types as HMM
import Biobase.SElab.Model.Types
-- | Filter a model after the header, not the body, has been parsed.
type PreFilterFun = Text -> Accession () -> Either (HMM ()) CM -> Bool
-- | Filter a model after the full model has been parsed. This is the same
-- type as @PreFilterFun@.
type PostFilterFun = Text -> Accession () -> Either (HMM ()) CM -> Bool
-- TODO use a builder?
newtype Log = Log { getLog :: Text }
deriving (Monoid,IsString)
-- | The type of logger we use
type Logger m = WriterT Log m
-- | Combine CMs with their HMMs. Assumes that each CM is followed by its
-- HMM.
attachHMMs :: (Monad m) => PP.Producer Model (Logger m) r -> PP.Producer CM (Logger m) ((), PP.Producer Model (Logger m) r)
attachHMMs = PP.parsed go where
go :: (Monad m) => PP.Parser Model (Logger m) (Either () CM)
go = do
mcm <- PP.draw
case mcm of
Nothing -> return $ Left ()
Just (Left hmm) -> do
lift . tell . Log $ "HMM: " <> (hmm^.HMM.name) <> " is an orphan\n"
go
Just (Right cm) -> do
mhm <- PP.draw
case mhm of
Nothing -> do
lift . tell . Log $ "CM: " <> (cm^.CM.name) <> " has no attached HMM and the stream is finished\n"
return $ Right cm
-- TODO actually check that these belong together
Just (Left hm) | (cm^.CM.name) == (hm^.HMM.name) -> do
-- cm and hmm belong together
return . Right $ set hmm (over HMM.accession retagAccession hm) cm
-- The HMM doesn't belong to our CM
Just (Left hm) -> do
lift . tell . Log $ "CM: " <> (cm^.CM.name) <> " and HMM: " <> (hm^.HMM.name) <> " do not belong together, dropping the HMM from the stream\n"
return . Right $ cm
Just (Right dup) -> do
lift . tell . Log $ "CM: " <> (cm^.CM.name) <> " has no attached HMM\n"
PP.unDraw $ Right dup
go
-- | Parses @HMM@ and @CM@ models from Rfam. The filtering function takes
-- the model name and accession and allows for premature termination of the
-- parsing of the current model.
parseSelectively :: (Monad m)
=> PreFilterFun
-- ^ filter function for premodels
-> PostFilterFun
-- ^ filter function for full models
-> PP.Producer ByteString (Logger m) r
-- -> PP.Producer (Maybe (Either (HMM xfam) CM)) (Logger m) (r, PP.Producer ByteString (Logger m) r)
-> PP.Producer Model (Logger m) ((), PP.Producer ByteString (Logger m) r)
parseSelectively preFltr postFltr p
= PP.parsed go p
P.>-> P.concat
P.>-> P.filter (\mdl -> postFltr (mdl^.modelName) (mdl^.modelAccession) mdl)
where
-- | Parse either a CM or a HMM ...
go = do
p <- zoom (splitKeepEnd "//\n") parseMdl
-- TODO can be simplified now
case p of
Left () -> return $ Left ()
Right x -> return $ Right x
-- parse models.
-- TODO @Either ()@ should become @Either (Maybe Error)@ and only @Left
-- Nothing@ will be error-free stop.
handleError err = do
da <- PP.drawAll
lift . tell . Log $ "could not parse:\n"
lift . tell . Log $ T.pack $ show err
lift . tell . Log . decodeUtf8 $ BS.concat da
lift . tell . Log $ "\n"
return $ Left ()
parseMdl :: Monad m => PP.StateT (PP.Producer ByteString (Logger m) x) (Logger m) (Either () (Maybe (Either (HMM ()) CM)))
parseMdl = do
-- if @pre@ is Nothing, the underlying producer is exhausted.
-- if @pre@ is @Just $ Left x@, then we have a parse error.
-- if @pre@ is @Just $ Right y@, then we have a successful parse. In
-- this case, @y@ is either a @Left hmm@ or a @Right cm@.
pre <- PA.parse $ (Left <$> parsePreHMM) <|> (Right <$> parsePreCM)
case pre of
-- we have nothing left to parse and indicate this now
Nothing -> return $ Left ()
Just (Left err) -> handleError err
Just (Right mdl) -> if preFltr (mdl^.modelName) (mdl^.modelAccession) mdl
then do
case mdl of
Left hmm -> do
h <- PA.parse $ parseHMMBody hmm
case h of
Nothing -> handleError "premature end of parsing in hmm body"
Just (Left err) -> handleError err
Just (Right hh) -> do
da <- PP.drawAll
-- TODO check if @da@ is empty?
return $ Right $ Just $ Left hh
Right cm -> do
c <- PA.parse $ parseCMBody cm
case c of
Nothing -> handleError "premature end of parsing in cm body"
Just (Left err) -> handleError err
Just (Right d) -> do
da <- PP.drawAll
return . Right . Just $ Right d
else PP.skipAll >> (return . Right $ Nothing)
-- | Keep all models
keepAllModels _ _ _ = True
-- | Load a number of models from file. Including pre- and full-model
-- filtering.
fromFile
:: PreFilterFun
-- ^ filter premodels before they are fully parsed. Full parsing is
-- costly. Use @\name acc hmmOrcm -> True@ if all models should be
-- loaded.
-> PostFilterFun
-- ^ Filter full models before they are combined into the CM-HMM pair.
-> Bool
-- ^ If true, than any error during parsing means termination of the
-- program.
-> FilePath
-- ^ input file name. Can be @-@ for stdin. If a file and the file ends
-- with @.gz@, the file is uncompressed on the ly.
-> IO [CM]
fromFile preFltr postFltr stopOnError fp
| fp == "-" = parse (PB.fromHandle stdin)
| takeExtension fp == ".gz" = withFile fp ReadMode $ \hdl -> parse (PG.decompress $ PB.fromHandle hdl)
| otherwise = withFile fp ReadMode $ \hdl -> parse (PB.fromHandle hdl)
where
parse source = do
((xs,((),rmdr)),log) <- runWriterT . P.toListM' $ attachHMMs $ parseSelectively preFltr postFltr source
-- TODO log should be empty
-- TODO rmdr should be empty
let Log l = log
unless (T.null l) $ do
T.putStrLn l
T.putStrLn "There have been errors parsing the models!"
when stopOnError $ do
error "stopping here!"
return xs