{-| Data type for a more general sequence location consiting of
potentially disjoint ranges of positions on the sequence.

Throughout, /sequence position/ refers to a 'Pos.Pos' which includes a
strand.  An index into a sequence is referred to as an /offset/, and
is generally of type 'Offset'.

 -}

module Bio.Location.Location ( 
  -- * Sequence locations
  Loc(..)

  -- * Locations and positions
  , bounds, length, startPos, endPos, posInto, posOutof, isWithin, overlaps

  -- * Extracting subsequences
  , seqData, seqDataPadded

  -- * Transforming locations
  , extend

  -- * Displaying locations
  , display
  ) where 

import Prelude hiding (length)

import Control.Arrow ((***))
import Control.Monad.Error
import qualified Data.ByteString.Lazy.Char8 as LBS
import Data.List (intercalate)

import qualified Bio.Location.ContigLocation as CLoc
import qualified Bio.Location.Position as Pos
import Bio.Location.Strand
import Bio.Sequence.SeqData

-- | General (disjoint) sequence region consisting of a concatenated
--   set of contiguous regions (see 'CLoc.ContigLoc').
newtype Loc = Loc [CLoc.ContigLoc] deriving (Eq, Ord, Show)

instance Stranded Loc where
    revCompl (Loc contigs) = Loc $ reverse $ map revCompl contigs

-- | Returns the length of the region
length :: Loc -> Offset
length (Loc contigs) = sum $ map CLoc.length contigs

-- | The bounds of a sequence location.  This is a pair consisting of
-- the lowest and highest sequence offsets covered by the region.  The
-- bounds ignore the strand of the sequence location, and the first
-- element of the pair will always be lower than the second.  Even if
-- the positions in the location do not run monotonically through the
-- location, the overall lowest and highest sequence offsets are returned.
bounds :: Loc -> (Offset, Offset)
bounds (Loc []) = error "locBounds on zero-contig Loc"
bounds (Loc contigs) = (minimum *** maximum) $ unzip $ map CLoc.bounds contigs

-- | Sequence position of the start of the location.  This is the 5'
-- end on the location strand, which will have a higher offset than
-- 'endPos' if the location is on the 'RevCompl' strand.
startPos :: Loc -> Pos.Pos
startPos (Loc [])      = error "startPos: zero-contig Loc"
startPos (Loc contigs) = CLoc.startPos $ head contigs

-- | Sequence position of the end of the location, as described in 'startPos'.
endPos :: Loc -> Pos.Pos
endPos (Loc [])      = error "endPos: zero-contig Loc"
endPos (Loc contigs) = CLoc.endPos $ last contigs

-- | Extract the nucleotide 'SeqData' for the sequence location.  If
-- any part of the location lies outside the bounds of the sequence,
-- an error results.
seqData :: (Error e, MonadError e m) => SeqData -> Loc -> m SeqData
seqData sequ (Loc contigs)
    = liftM LBS.concat $ mapM (CLoc.seqData sequ) contigs

-- | As 'seqData', extract the nucleotide subsequence for the
-- location.  Any positions in the location lying outside the bounds
-- of the sequence are returned as @N@ rather than producing an error.
seqDataPadded :: SeqData -> Loc -> SeqData
seqDataPadded sequ (Loc contigs)
    = LBS.concat $ map (CLoc.seqDataPadded sequ) contigs

-- | Given a sequence position and a sequence location relative to the
-- same sequence, compute a new position representing the original
-- position relative to the subsequence defined by the location.  If
-- the sequence position lies outside of the sequence location,
-- @Nothing@ is returned; thus, the offset of the new position will
-- always be in the range @[0, length cloc - 1]@.
--
-- When the sequence positions in the location are not monotonic,
-- there may be multiple possible posInto solutions.  That is, if the
-- same outer sequence position is covered by two different contiguous
-- blocks of the location, then it would have two possible sequence
-- positions relative to the location. In this case, the position
-- 5\'-most in the location orientation is returned.
posInto :: Pos.Pos -> Loc -> Maybe Pos.Pos
posInto seqpos (Loc contigs) = posIntoContigs seqpos contigs

posIntoContigs :: Pos.Pos -> [CLoc.ContigLoc] -> Maybe Pos.Pos
posIntoContigs _ [] = Nothing
posIntoContigs seqpos (contig@(CLoc.ContigLoc _ len _):rest)
    = case CLoc.posInto seqpos contig of
        just@(Just _) -> just
        Nothing -> liftM (flip Pos.slide len) $ posIntoContigs seqpos rest

-- | Given a sequence location and a sequence position within that
-- location, compute a new position representing the original position
-- relative to the outer sequence.  If the sequence position lies
-- outside the location, @Nothing@ is returned.  
-- 
-- This function inverts 'posInto' when the sequence position lies
-- within the position is actually within the location.  Due to the
-- possibility of redundant location-relative positions for a given
-- absolute position, 'posInto' does not necessary invert 'posOutof'
posOutof :: Pos.Pos -> Loc -> Maybe Pos.Pos
posOutof pos (Loc contigs) = posOutofContigs pos contigs

posOutofContigs :: Pos.Pos -> [CLoc.ContigLoc] -> Maybe Pos.Pos
posOutofContigs _ [] = Nothing
posOutofContigs seqpos (contig@(CLoc.ContigLoc _ len _):rest)
    = case CLoc.posOutof seqpos contig of
        just@(Just _) -> just
        Nothing -> posOutofContigs (Pos.slide seqpos $ negate len) rest

-- | Returns a sequence location produced by extending the original
-- location on each end, based on a pair of (/5\' extension/, /3\'
-- extension/).  These add contiguous positions to the 5\' and 3\'
-- ends of the original location.  The 5\' extension is applied to the
-- 5\' end of the location on the location strand; if the location is
-- on the 'RevCompl' strand, the 5\' end will have a higher offset
-- than the 3\' end and this offset will increase by the amount of the
-- 5\' extension.  Similarly, the 3\' extension is applied to the 3\'
-- end of the location.
extend :: (Offset, Offset) -- ^ (5' extension, 3' extension)
       -> Loc -> Loc
extend _ (Loc []) = error "extendLoc on zero-contig Loc"
extend (ext5, ext3) (Loc contigs) = Loc $ case extendContigs3 contigs of
                                            [] -> error "Empty contigs after extendContigs3"
                                            (cfirst:crest) -> (CLoc.extend (ext5, 0) cfirst):crest
    where extendContigs3 [] = error "Empty contigs in extendContigs3"
          extendContigs3 [clast] = [CLoc.extend (0, ext3) clast]
          extendContigs3 (contig:crest) = contig : extendContigs3 crest

-- | Returns @True@ when  a sequence position lies within a sequence
-- location on the same sequence, and occupies the same strand.
isWithin :: Pos.Pos -> Loc -> Bool
isWithin seqpos (Loc contigs) = or $ map (CLoc.isWithin seqpos) contigs

overlappingContigs :: Loc -> Loc -> [(CLoc.ContigLoc, CLoc.ContigLoc)]
overlappingContigs (Loc contigs1) (Loc contigs2) 
    = filter (uncurry CLoc.overlaps) [(c1, c2) | c1 <- contigs1, c2 <- contigs2 ]

-- | Returns @True@ when two sequence locations overlap at any
-- position.
overlaps :: Loc -> Loc -> Bool
overlaps l1 l2 = not $ null $ overlappingContigs l1 l2

-- | Display a human-friendly, zero-based representation of a sequence location.
display :: Loc -> String
display (Loc contigs) = intercalate ";" $ map CLoc.display contigs