samtools-iteratee-0.1: src/BamFilter.hs

module Main
       where

import Control.Applicative
import Control.Monad.Reader
import Data.Maybe

import System.Console.GetOpt
import System.Environment
import System.IO

import qualified Data.Iteratee as Iter

import qualified Bio.SamTools.Bam as Bam
import qualified Bio.SamTools.Iteratee as Bam

main :: IO ()
main = getArgs >>= handleOpt . getOpt RequireOrder optDescrs
    where handleOpt (_,    _,         errs@(_:_)) = usage (unlines errs)
          handleOpt (args, [bam], []) = either usage (doBamFilter bam) $ argsToConf args
          handleOpt (_,    _,     []) = usage "Specify exactly one BAM input"
          usage errs = do prog <- getProgName
                          hPutStr stderr $ usageInfo prog optDescrs
                          hPutStrLn stderr errs

doBamFilter :: FilePath -> Conf -> IO ()
doBamFilter bam conf = Bam.withBamInFile bam $ \hin ->
  Bam.withBamOutFile (confOutput conf)  (Bam.inHeader hin) $ \hout ->
  let bamiter = Iter.joinI $ Iter.filter (wanted conf) $ Iter.mapM_ (Bam.put1 hout)
  in Bam.enumInHandle hin bamiter >>= Iter.run
  
wanted :: Conf -> Bam.Bam1 -> Bool
wanted conf | confPerfect conf = maybe False (== 0) . Bam.nMismatch
            | otherwise = const True
                                    
data Conf = Conf { confOutput :: !FilePath
                 , confPerfect :: !Bool
                 } deriving (Show)

data Arg = ArgOutput { unArgOutput :: !String }
         | ArgPerfect
         deriving (Show, Read, Eq, Ord)

argOutput :: Arg -> Maybe String
argOutput (ArgOutput del) = Just del
argOutput _ = Nothing

optDescrs :: [OptDescr Arg]
optDescrs = [ Option ['o'] ["output"]  (ReqArg ArgOutput "OUTFILE") "Output filename"
            , Option []    ["perfect"] (NoArg ArgPerfect)           "Filter perfect (NM:i:0) reads"
            ]

argsToConf :: [Arg] -> Either String Conf
argsToConf = runReaderT conf
    where conf = Conf <$> 
                 findOutput <*>
                 (ReaderT $ return . elem ArgPerfect)
          findOutput = ReaderT $ maybe (Left "No output file") return . listToMaybe . mapMaybe argOutput