cobot-io 0.1.2.10 → 0.1.3.0

raw patch · 17 files changed

+1034/−550 lines, 17 filesPVP: major bump suggested

API removals or changes: PVP suggests a major version bump

API changes (from Hackage documentation)

Files

• ChangeLog.md +11/−0
• cobot-io.cabal +5/−2
• src/Bio/MAE.hs +22/−8
• src/Bio/MAE/Type.hs +9/−0
• src/Bio/MMTF.hs +11/−8
• src/Bio/PDB.hs +63/−42
• src/Bio/PDB/BondRestoring.hs +196/−0
• src/Bio/PDB/Functions.hs +20/−0
• src/Bio/PDB/Parser.hs +8/−3
• src/Bio/PDB/Reader.hs +10/−9
• src/Bio/PDB/Type.hs +3/−0
• src/Bio/Structure.hs +18/−12
• test/MAESpec.hs +85/−64
• test/MMTFSpec.hs +4/−4
• test/PDBParserSpec.hs +446/−0
• test/PDBSpec.hs +119/−398
• test/Spec.hs +4/−0

@@ -2,6 +2,17 @@  ## [Unreleased] +## [0.1.3.0] - 2020-03-27+### Added+- Residue index in `Structure`.+- Atom input index in `Structure`.+- Bond restoring for PDB.+- Tests for PDB -> Model conversion.+### Changed+- GlobalID now 0-based in mae, PDB, and MMTF.+### Fixed+- A lot of things.+ ## [0.1.2.10] - 2020-03-27 ### Added - Lenses for `Structure`.

@@ -4,10 +4,10 @@ -- -- see: https://github.com/sol/hpack ----- hash: 9e1717465825e1ea9fa611d37fb048cda60ddcea37b8afd16ba77a858416e927+-- hash: a5a2fbb335e61a22b55a3d74cd08557730ea19a8b93005c693a4576dad82e498  name:           cobot-io-version:        0.1.2.10+version:        0.1.3.0 synopsis:       Biological data file formats and IO description:    Please see the README on GitHub at <https://github.com/less-wrong/cobot-io#readme> category:       Bio@@ -50,6 +50,8 @@       Bio.MMTF.MessagePack       Bio.MMTF.Type       Bio.PDB+      Bio.PDB.BondRestoring+      Bio.PDB.Functions       Bio.PDB.Parser       Bio.PDB.Reader       Bio.PDB.Type@@ -103,6 +105,7 @@       MAEParserSpec       MAESpec       MMTFSpec+      PDBParserSpec       PDBSpec       SequenceSpec       StructureSpec

@@ -1,5 +1,5 @@-{-# LANGUAGE AllowAmbiguousTypes #-}-{-# LANGUAGE CPP                 #-}+{-# LANGUAGE AllowAmbiguousTypes  #-}+{-# LANGUAGE CPP                  #-} {-# OPTIONS_GHC -fno-warn-orphans #-}  module Bio.MAE@@ -8,13 +8,15 @@   , Table (..)   , fromFile   , fromText+  , modelsFromMaeText+  , modelsFromMaeFile   , maeP   ) where  import           Bio.MAE.Parser import           Bio.MAE.Type           (Block (..), FromMaeValue (..),                                          Mae (..), MaeValue (..), Table (..))-import           Bio.Structure          (Atom (..), Bond (..), Chain (..),+import           Bio.Structure          (Atom (..), Bond (..), Chain (..), Model (..),                                          GlobalID (..), LocalID (..),                                          Model (..), Residue (..),                                          SecondaryStructure (..),@@ -50,6 +52,12 @@ fromText :: Text -> Either Text Mae fromText = first T.pack . parseOnly maeP +modelsFromMaeFile :: (MonadIO m) => FilePath -> m (Either Text (Vector Model))+modelsFromMaeFile = liftIO . fmap modelsFromMaeText . TIO.readFile++modelsFromMaeText :: Text -> Either Text (Vector Model)+modelsFromMaeText maeText = modelsOf <$> fromText maeText+ instance StructureModels Mae where   modelsOf Mae{..} = V.fromList $ fmap blockToModel blocks     where@@ -117,18 +125,23 @@                   groupedByResidues = toGroupsOn by group                   residues          = V.fromList $ fmap groupToResidue groupedByResidues -                  by :: Int -> (Int, Text)-                  by i = (unsafeGetFromContents "i_m_residue_number" i, getFromContents defaultChainName "s_m_insertion_code" i)+                  by :: Int -> (Int, Char)+                  by i = (unsafeGetFromContents "i_m_residue_number" i, getFromContents defaultInsertionCode "s_m_insertion_code" i)                defaultChainName :: Text               defaultChainName = "A" +              defaultInsertionCode :: Char+              defaultInsertionCode = ' '+               groupToResidue :: [Int] -> Residue               groupToResidue []            = error "Group that is result of List.groupBy can't be empty."-              groupToResidue group@(h : _) = Residue name atoms (V.fromList localBonds) secondary chemCompType+              groupToResidue group@(h : _) = Residue name residueNumber insertionCode atoms (V.fromList localBonds) secondary chemCompType                 where-                  name  = stripQuotes $ unsafeGetFromContents "s_m_pdb_residue_name" h-                  atoms = V.fromList $ fmap indexToAtom group+                  name          = stripQuotes $ unsafeGetFromContents "s_m_pdb_residue_name" h+                  residueNumber = unsafeGetFromContents "i_m_residue_number" h+                  insertionCode = unsafeGetFromContents "s_m_insertion_code" h+                  atoms         = V.fromList $ fmap indexToAtom group                    localInds     = [0 .. length group - 1]                   globalToLocal = M.fromList $ zip group localInds@@ -146,6 +159,7 @@                indexToAtom :: Int -> Atom               indexToAtom i = Atom (GlobalID i)+                                   (i + 1)                                    (stripQuotes $ getFromContentsI "s_m_pdb_atom_name")                                    (elIndToElement M.! getFromContentsI "i_m_atomic_number")                                    coords

@@ -9,6 +9,7 @@ import           Data.Map.Strict (Map) import           Data.Maybe      (fromJust) import           Data.Text       (Text)+import qualified Data.Text       as T (head, null, dropAround)  data Mae = Mae { version :: Text                , blocks  :: [Block]@@ -58,3 +59,11 @@     fromMaeValue :: MaeValue -> Maybe Text     fromMaeValue (StringMaeValue t) = Just t     fromMaeValue _                  = Nothing++instance FromMaeValue Char where+    fromMaeValue :: MaeValue -> Maybe Char+    fromMaeValue (StringMaeValue t) = Just $ if T.null t then ' ' else T.head $ stripQuotes t+    fromMaeValue _                  = Nothing+         +stripQuotes :: Text -> Text+stripQuotes = T.dropAround (== '"')

@@ -41,6 +41,7 @@                  decode (getResponseBody resp)  instance StructureModels MMTF where+    -- TODO: add global bonds     modelsOf m = l2v (flip Model empty . l2v <$> zipWith (zipWith Chain) chainNames chainResis)       where         chainsCnts = fromIntegral <$> toList (chainsPerModel (model m))@@ -70,7 +71,8 @@                              in  (end, mkAtom <$> zip4 cl nl el ics)          mkResidue :: (GroupType, SecondaryStructure, [Atom]) -> Residue-        mkResidue (gt, ss, atoms') = Residue (gtGroupName gt) (l2v atoms')+        -- TODO: support residue number here+        mkResidue (gt, ss, atoms') = Residue (gtGroupName gt) (-1) ' ' (l2v atoms')                                              (mkBonds (gtBondAtomList gt) (gtBondOrderList gt))                                               ss (gtChemCompType gt) @@ -87,13 +89,14 @@                                      z = zCoordList (atom m)                                      o = occupancyList (atom m)                                      b = bFactorList (atom m)-                                 in  Atom (GlobalID $ fromIntegral (i ! idx))-                                           n-                                           e-                                           (V3 (x ! idx) (y ! idx) (z ! idx))-                                           fc-                                           (b ! idx)-                                           (o ! idx)+                                 in  Atom (GlobalID idx)+                                          (fromIntegral $ i ! idx)+                                          n+                                          e+                                          (V3 (x ! idx) (y ! idx) (z ! idx))+                                          fc+                                          (b ! idx)+                                          (o ! idx)          cutter :: [Int] -> [a] -> [[a]]         cutter []     []    = []

@@ -1,61 +1,82 @@ {-# OPTIONS_GHC -fno-warn-orphans #-} module Bio.PDB-  (+  ( modelsFromPDBText+  , modelsFromPDBFile   ) where -import qualified Bio.PDB.Type  as PDB+import qualified Bio.PDB.Type           as PDB+import           Bio.PDB.Reader         (fromTextPDB, PDBWarnings)+import           Bio.PDB.BondRestoring  (restoreModelGlobalBonds, restoreChainLocalBonds, residueID)+import           Bio.PDB.Functions      (groupChainByResidue) import           Bio.Structure -import           Control.Arrow ((&&&))-import           Data.Coerce   (coerce)-import           Data.Foldable (Foldable (..))-import           Data.Text     as T (Text, singleton, unpack)-import qualified Data.Vector   as V-import           Linear.V3     (V3 (..))+import           Control.Arrow          ((&&&))+import           Control.Monad.IO.Class (MonadIO, liftIO) +import           Data.Text              as T (Text, singleton, unpack, strip)+import           Data.Text.IO           as TIO (readFile)+import           Data.Map               (Map)+import qualified Data.Map               as M ((!), fromList)+import qualified Data.Vector            as V+import           Data.List              (sort)+import           Data.Maybe             (fromMaybe)++import           Text.Read              (readMaybe)++import           Linear.V3              (V3 (..))+ instance StructureModels PDB.PDB where     modelsOf PDB.PDB {..} = fmap mkModel models       where         mkModel :: PDB.Model -> Model-        mkModel = flip Model V.empty . fmap mkChain--        mkChain :: PDB.Chain -> Chain-        mkChain = uncurry Chain . (mkChainName &&& mkChainResidues)+        mkModel model = Model (fmap mkChain model) (restoreModelGlobalBonds atomSerialToNilBasedIndex model)+          where+            atomSerialToNilBasedIndex :: Map Int Int+            atomSerialToNilBasedIndex = M.fromList $ allModelAtomSerials `zip` [0..] -        mkChainName :: PDB.Chain -> Text-        mkChainName = T.singleton . PDB.atomChainID . safeFirstAtom+            allModelAtomSerials :: [Int]+            allModelAtomSerials = sort . V.toList . fmap PDB.atomSerial . V.concat $ V.toList model -        mkChainResidues :: PDB.Chain -> V.Vector Residue-        mkChainResidues = V.fromList . fmap mkResidue . flip groupByResidue [] . pure . toList+            mkChain :: PDB.Chain -> Chain+            mkChain = uncurry Chain . (mkChainName &&& mkChainResidues) -        -- can be rewritten with sortOn and groupBy-        groupByResidue :: [[PDB.Atom]] -> [PDB.Atom] -> [[PDB.Atom]]-        groupByResidue res []       = res-        groupByResidue [] (x : xs)  = groupByResidue [[x]] xs-        groupByResidue res@(lastList : resultTail) (x : xs)-          | (PDB.atomResSeq x, PDB.atomICode x) == (PDB.atomResSeq (head lastList), PDB.atomICode (head lastList))-                                              = groupByResidue ((x : lastList) : resultTail) xs-          | otherwise                         = groupByResidue ([x] : res) xs+            mkChainName :: PDB.Chain -> Text+            mkChainName = T.singleton . PDB.atomChainID . safeFirstAtom -        safeFirstAtom :: V.Vector PDB.Atom -> PDB.Atom-        safeFirstAtom arr | V.length arr > 0 = arr V.! 0-                          | otherwise        = error "Could not pick first atom"+            mkChainResidues :: PDB.Chain -> V.Vector Residue+            mkChainResidues chain = V.fromList . fmap (mkResidue (restoreChainLocalBonds chain)) $ groupChainByResidue chain +            safeFirstAtom :: V.Vector PDB.Atom -> PDB.Atom+            safeFirstAtom arr | V.length arr > 0 = arr V.! 0+                              | otherwise        = error "Could not pick first atom"+            +            mkResidue :: Map Text (V.Vector (Bond LocalID)) -> [PDB.Atom] -> Residue+            mkResidue _ []    = error "Cound not make residue from empty list"+            mkResidue localBondsMap atoms' = Residue (T.strip $ PDB.atomResName firstResidueAtom)+                                                     (PDB.atomResSeq firstResidueAtom)+                                                     (PDB.atomICode firstResidueAtom)+                                                     (V.fromList $ mkAtom <$> atoms')+                                                     (localBondsMap M.! residueID firstResidueAtom)+                                                     Undefined -- now we do not read secondary structure+                                                     ""        -- chemical component type?!+              where+                firstResidueAtom = head atoms' -        mkResidue :: [PDB.Atom] -> Residue-        mkResidue []     = error "Cound not make residue from empty list"-        mkResidue atoms' = Residue (PDB.atomResName . head $ atoms')-                                   (V.fromList $ mkAtom <$> atoms')-                                   V.empty   -- now we do not read bonds-                                   Undefined -- now we do not read secondary structure-                                   ""        -- chemical component type?!+            mkAtom :: PDB.Atom -> Atom+            mkAtom PDB.Atom{..} = Atom (GlobalID $ atomSerialToNilBasedIndex M.! atomSerial)+                                       atomSerial+                                       (T.strip atomName)+                                       atomElement+                                       (V3 atomX atomY atomZ)+                                       (fromMaybe 0 . readMaybe $ T.unpack atomCharge)+                                       atomTempFactor+                                       atomOccupancy +modelsFromPDBFile :: (MonadIO m) => FilePath -> m (Either Text ([PDBWarnings], V.Vector Model))+modelsFromPDBFile = liftIO . fmap modelsFromPDBText . TIO.readFile -        mkAtom :: PDB.Atom -> Atom-        mkAtom PDB.Atom{..} = Atom (coerce atomSerial)-                                   atomName-                                   atomElement-                                   (V3 atomX atomY atomZ)-                                   (read $ T.unpack atomCharge)-                                   atomTempFactor-                                   atomOccupancy+modelsFromPDBText :: Text -> Either Text ([PDBWarnings], V.Vector Model)+modelsFromPDBText pdbText = do+  (warnings, parsedPDB) <- fromTextPDB pdbText+  let models = modelsOf parsedPDB+  pure (warnings, models)

@@ -0,0 +1,196 @@+{-# LANGUAGE OverloadedStrings #-}+{-# LANGUAGE TupleSections #-}++module Bio.PDB.BondRestoring+  ( restoreModelGlobalBonds+  , restoreModelLocalBonds+  , restoreChainLocalBonds+  , residueID+  ) where++import qualified Bio.PDB.Type as PDB  (Atom(..))+import           Bio.PDB.Functions    (groupChainByResidue)+import           Bio.Structure        (Bond (..), GlobalID (..), LocalID (..))++import           Data.Vector          (Vector)+import qualified Data.Vector as V     (fromList, toList)+import           Data.List            (find, sort)+import           Data.Text            (Text)+import qualified Data.Text as T       (strip, pack, unpack)+import           Data.Map.Strict      (Map, (!?), (!))+import qualified Data.Map.Strict as M (fromList)+import           Data.Maybe           (catMaybes)++import           Linear.Metric        (distance)+import           Linear.V3            (V3(..))++import           Control.Monad        (guard)++residueID :: PDB.Atom -> Text+residueID PDB.Atom{..} = T.pack (show atomChainID) <> T.pack (show atomResSeq) <> T.pack (show atomICode)++restoreModelLocalBonds :: Vector (Vector PDB.Atom) -> Map Text (Vector (Bond LocalID))+restoreModelLocalBonds = M.fromList . concatMap restoreChainLocalBonds'++restoreChainLocalBonds :: Vector PDB.Atom -> Map Text (Vector (Bond LocalID))+restoreChainLocalBonds = M.fromList . restoreChainLocalBonds'++restoreChainLocalBonds' :: Vector PDB.Atom -> [(Text, Vector (Bond LocalID))]+restoreChainLocalBonds' chainAtoms = residueIDToLocalBonds+  where+    residueIDToLocalBonds :: [(Text, Vector (Bond LocalID))]+    residueIDToLocalBonds = do+      (residueAtoms, residueBonds) <- zip chainAtomsGroupedByResidue intraResidueGlobalBonds+      let localBonds = V.fromList $ convertGlobalsToLocals residueAtoms residueBonds+      let _residueID = residueID $ head residueAtoms+      pure (_residueID, localBonds)+    +    intraResidueGlobalBonds :: [[Bond GlobalID]]+    intraResidueGlobalBonds = fmap restoreIntraResidueBonds chainAtomsGroupedByResidue+    +    chainAtomsGroupedByResidue :: [[PDB.Atom]]+    chainAtomsGroupedByResidue = groupChainByResidue chainAtoms+    +    convertGlobalsToLocals :: [PDB.Atom] -> [Bond GlobalID] -> [Bond LocalID]+    convertGlobalsToLocals residueAtoms = map convertGlobalToLocal+      where+        convertGlobalToLocal :: Bond GlobalID -> Bond LocalID+        convertGlobalToLocal (Bond (GlobalID from) (GlobalID to) order) = +          Bond (LocalID $ globalToLocalIdxMap ! from) (LocalID $ globalToLocalIdxMap ! to) order+        +        globalToLocalIdxMap :: Map Int Int+        globalToLocalIdxMap = M.fromList $ zip sortedGlobalIndices [0..]+        +        sortedGlobalIndices :: [Int]+        sortedGlobalIndices = map PDB.atomSerial $ sort residueAtoms+++restoreModelGlobalBonds :: Map Int Int -> Vector (Vector PDB.Atom) -> Vector (Bond GlobalID)+restoreModelGlobalBonds atomSerialToNilBasedIndex chains = convertGlobalIDs atomSerialToNilBasedIndex . V.fromList $ _intraResidueBonds ++ peptideBonds ++ disulfideBonds+  where+    convertGlobalIDs :: Map Int Int -> Vector (Bond GlobalID) -> Vector (Bond GlobalID)+    convertGlobalIDs mapping = reindexBonds (\(GlobalID v) -> GlobalID $ mapping ! v)+    +    reindexBonds :: (a -> a) -> Vector (Bond a) -> Vector (Bond a)+    reindexBonds convertID = fmap (\(Bond from to order) -> Bond (convertID from) (convertID to) order)++    chainAtomsGroupedByResidue :: Vector [[PDB.Atom]]+    chainAtomsGroupedByResidue = fmap groupChainByResidue chains+    +    _intraResidueBonds :: [Bond GlobalID]+    _intraResidueBonds = concatMap restoreChainIntraResidueBonds chainAtomsGroupedByResidue+    +    peptideBonds :: [Bond GlobalID]+    peptideBonds = concatMap restoreChainPeptideBonds chainAtomsGroupedByResidue+    +    disulfideBonds :: [Bond GlobalID]+    disulfideBonds = restoreDisulfideBonds . concat $ V.toList chainAtomsGroupedByResidue++restoreDisulfideBonds :: [[PDB.Atom]] -> [Bond GlobalID]+restoreDisulfideBonds atomsGroupedByResidue = do+  atom1 <- cystineSulfur+  atom2 <- cystineSulfur+  guard (PDB.atomSerial atom1 < PDB.atomSerial atom2)+  guard $ distance (coords atom1) (coords atom2) < sulfidicBondMaxLength+  pure $ Bond (GlobalID $ PDB.atomSerial atom1) (GlobalID $ PDB.atomSerial atom2) 1+  where+    cystineSulfur :: [PDB.Atom]+    cystineSulfur = filter (("SG" ==) . T.strip . PDB.atomName) $ concat cystines+    cystines :: [[PDB.Atom]]+    cystines = filter cystinePredicate atomsGroupedByResidue+    cystinePredicate :: [PDB.Atom] -> Bool+    cystinePredicate residue = any (("SG" ==) . T.strip . PDB.atomName) residue && all (("HG" /=) . T.strip . PDB.atomName) residue+    coords :: PDB.Atom -> V3 Float+    coords PDB.Atom{..} = V3 atomX atomY atomZ++sulfidicBondMaxLength :: Float+sulfidicBondMaxLength = 2.56++restoreChainPeptideBonds :: [[PDB.Atom]] -> [Bond GlobalID]+restoreChainPeptideBonds atomsGroupedByResidue = restoreChainPeptideBonds' atomsGroupedByResidue []+  where+    restoreChainPeptideBonds' :: [[PDB.Atom]] -> [Bond GlobalID] -> [Bond GlobalID]+    restoreChainPeptideBonds' [] acc = acc+    restoreChainPeptideBonds' [_] acc = acc+    restoreChainPeptideBonds' (residue1:residue2:residues) acc = +      restoreChainPeptideBonds' (residue2:residues) (constructBond residue1 residue2 : acc)++    constructBond :: [PDB.Atom] -> [PDB.Atom] -> Bond GlobalID+    constructBond residue1 residue2 = Bond (GlobalID $ getAtomIndex residue1 "C") (GlobalID $ getAtomIndex residue2 "N") 1+    +    getAtomIndex :: [PDB.Atom] -> Text -> Int+    getAtomIndex atoms atomNameToFind = case find ((atomNameToFind ==) . T.strip . PDB.atomName) atoms of +      Just PDB.Atom{..} -> atomSerial+      Nothing           -> error ("Atom with name " ++ T.unpack atomNameToFind ++ " wasn't found in residue " ++ residueId atoms ++ ", chain: " ++ chainId atoms)+    residueId :: [PDB.Atom] -> String+    residueId [] = error "cobot-io: it's impossible to form a residue ID on a residue with no atoms"+    residueId (PDB.Atom{..}:_) = T.unpack atomResName ++ show atomResSeq ++ show atomICode+    +    chainId :: [PDB.Atom] -> String+    chainId [] = error "cobot-io: it's impossible to get a chain ID on a chain with no atoms"+    chainId (PDB.Atom{..}:_) = show atomChainID+    ++restoreChainIntraResidueBonds :: [[PDB.Atom]] -> [Bond GlobalID]+restoreChainIntraResidueBonds = concatMap restoreIntraResidueBonds++restoreIntraResidueBonds :: [PDB.Atom] -> [Bond GlobalID]+restoreIntraResidueBonds residueAtoms = catMaybes $ constructBond <$> residueBonds+  where+    -- TODO: support bond order somehow+    constructBond :: (Text, Text) -> Maybe (Bond GlobalID)+    constructBond (fromAtomName, toAtomName) = Bond <$> constructGlobalID fromAtomName <*> constructGlobalID  toAtomName <*> Just 1+    +    constructGlobalID :: Text -> Maybe GlobalID+    constructGlobalID atomName = GlobalID <$> atomNameToIndex !? atomName+    +    atomNameToIndex :: Map Text Int+    atomNameToIndex = M.fromList $ (\PDB.Atom{..} -> (T.strip atomName, atomSerial)) <$> residueAtoms+    +    residueBonds :: [(Text, Text)]+    residueBonds = intraResidueBonds . T.strip . PDB.atomResName $ head residueAtoms++intraResidueBonds :: Text -> [(Text, Text)]+intraResidueBonds "NMA" = [("CA", "N")]+intraResidueBonds "ACE" = [("C", "O"), ("C", "CH3")]+intraResidueBonds residueName = backboneBonds ++ caCbBonds residueName ++ sideChainBonds residueName++backboneBonds :: [(Text, Text)]+backboneBonds = [("N", "CA"), ("CA", "C"), ("C", "O"), ("N", "H")] ++ [("C","OXT"), ("C","HXT")] ++ bwhMany [("N", ["H1", "H2", "H3"])]++caCbBonds :: Text -> [(Text, Text)]+caCbBonds aminoacid = case aminoacid of+  "GLY" -> bwhMany [("CA", ["HA2", "HA3"])]+  _     -> [("CA", "CB"), ("CA", "HA")]++sideChainBonds :: Text -> [(Text, Text)]+sideChainBonds "ALA" = bwhMany [("CB", ["HB1", "HB2", "HB3"])]+sideChainBonds "ARG" = [("CB", "CG"), ("CG", "CD"), ("CD", "NE"), ("NE", "CZ"), ("CZ", "NH2"), ("CZ", "NH1")] ++ bwhMany[("CB", ["HB3", "HB2"]), ("CG", ["HG3", "HG2"]), ("CD", ["HD3", "HD2"]), ("NE", ["HE"]), ("NH1", ["HH12", "HH11"]), ("NH2", ["HH22", "HH21"])]+sideChainBonds "ASN" = [("CB", "CG"), ("CG", "OD1"), ("CG", "ND2")] ++ bwhMany [("CB", ["HB3", "HB2"]), ("ND2", ["HD22", "HD21"])]+sideChainBonds "ASP" = [("CB", "CG"), ("CG", "OD1"), ("CG", "OD2")] ++ bwhMany [("CB", ["HB3", "HB2"]), ("OD2", ["HD2"])]+sideChainBonds "CYS" = [("CB", "SG")] ++ bwhMany [("CB", ["HB3", "HB2"]), ("SG", ["HG"])]+sideChainBonds "GLN" = [("CB", "CG"), ("CG", "CD"), ("CD", "OE1"), ("CD", "NE2")] ++ bwhMany [("CB", ["HB3", "HB2"]), ("CG", ["HG3", "HG2"]), ("NE2", ["HE22", "HE21"])]+sideChainBonds "GLU" = [("CB", "CG"), ("CG", "CD"), ("CD", "OE1"), ("CD", "OE2")] ++ bwhMany [("CB", ["HB3", "HB2"]), ("CG", ["HG3", "HG2"]), ("OE2", ["HE2"])]+sideChainBonds "GLY" = [] -- nothing+sideChainBonds "HIS" = [("CB", "CG"), ("CG", "ND1"), ("ND1", "CE1"), ("CE1", "NE2"), ("NE2", "CD2"), ("CD2", "CG")] ++ bwhMany [("CB", ["HB3", "HB2"]), ("ND1", ["HD1"]), ("CE1", ["HE1"]), ("NE2", ["HE2"]), ("CD2", ["HD2"])]+sideChainBonds "ILE" = [("CB", "CG1"), ("CB", "CG2"), ("CG1", "CD1")] ++ bwhMany [("CB", ["HB"]), ("CG1", ["HG13", "HG12"]), ("CG2", ["HG21", "HG22", "HG23"]), ("CD1", ["HD11", "HD12", "HD13"])]+sideChainBonds "LEU" = [("CB", "CG"), ("CG", "CD1"), ("CG", "CD2")] ++ bwhMany [("CB", ["HB3", "HB2"]), ("CG", ["HG"]), ("CD1", ["HD11", "HD12", "HD13"]), ("CD2", ["HD21", "HD22", "HD23"])]+sideChainBonds "LYS" = [("CB", "CG"), ("CG", "CD"), ("CD", "CE"), ("CE", "NZ")] ++ bwhMany [("CB", ["HB3", "HB2"]), ("CG", ["HG3", "HG2"]), ("CD", ["HD3", "HD2"]), ("CE", ["HE3", "HE2"]), ("NZ", ["HZ1", "HZ2", "HZ3"])]+sideChainBonds "MET" = [("CB", "CG"), ("CG", "SD"), ("SD", "CE")] ++ bwhMany [("CB", ["HB3", "HB2"]), ("CG", ["HG3", "HG2"]), ("CE", ["HE1", "HE2", "HE3"])]+sideChainBonds "PHE" = [("CB", "CG"), ("CG", "CD1"), ("CD1", "CE1"), ("CE1", "CZ"), ("CZ", "CE2"), ("CE2", "CD2"), ("CD2", "CG")] ++ bwhMany [("CB", ["HB3", "HB2"]), ("CD1", ["HD1"]), ("CE1", ["HE1"]), ("CZ", ["HZ"]), ("CE2", ["HE2"]), ("CD2", ["HD2"])]+sideChainBonds "PRO" = [("CB", "CG"), ("CG", "CD"), ("CD", "N")] ++ bwhMany [("CB", ["HB3", "HB2"]), ("CG", ["HG3", "HG2"]), ("CD", ["HD2", "HD3"])]+sideChainBonds "SER" = [("CB", "OG")] ++ bwhMany [("CB", ["HB3", "HB2"]), ("OG", ["HG"])]+sideChainBonds "THR" = [("CB", "OG1"), ("CB", "CG2")] ++ bwhMany [("CB", ["HB"]), ("OG1", ["HG1"]), ("CG2", ["HG21", "HG22", "HG23"])]+sideChainBonds "TRP" = [("CB", "CG"), ("CG", "CD1"), ("CD1", "NE1"), ("NE1", "CE2"), ("CE2", "CD2"), ("CD2", "CG"), ("CD2", "CE3"), ("CE3", "CZ3"), ("CZ3", "CH2"), ("CH2", "CZ2"), ("CZ2", "CE2")] ++ bwhMany [("CB", ["HB3", "HB2"]), ("CD1", ["HD1"]), ("NE1", ["HE1"]), ("CE3", ["HE3"]), ("CZ3", ["HZ3"]), ("CH2", ["HH2"]), ("CZ2", ["HZ2"])]+sideChainBonds "TYR" = [("CB", "CG"), ("CG", "CD1"), ("CD1", "CE1"), ("CE1", "CZ"), ("CZ", "CE2"), ("CE2", "CD2"), ("CD2", "CG"), ("CZ", "OH")] ++ bwhMany [("CB", ["HB3", "HB2"]), ("CD1", ["HD1"]), ("CE1", ["HE1"]), ("CE2", ["HE2"]), ("CD2", ["HD2"]), ("OH", ["HH"])]+sideChainBonds "VAL" = [("CB", "CG1"), ("CB", "CG2")] ++ bwhMany [("CB", ["HB"]), ("CG1", ["HG11", "HG12", "HG13"]), ("CG2", ["HG21", "HG22", "HG23"])]+sideChainBonds unknownResidue = error . T.unpack $ "cobot-io: we don't know what to do with residue " <> unknownResidue++bwhMany :: [(Text, [Text])] -> [(Text, Text)]+bwhMany = concatMap bwh++bwh :: (Text, [Text]) -> [(Text, Text)]+bwh = heavyAtomBondsWithHydrogens++heavyAtomBondsWithHydrogens :: (Text, [Text]) -> [(Text, Text)]+heavyAtomBondsWithHydrogens (heavyAtomName, hydrogenNames) = (heavyAtomName,) <$> hydrogenNames

@@ -0,0 +1,20 @@+module Bio.PDB.Functions+  ( groupChainByResidue+  ) where++import qualified Bio.PDB.Type as PDB (Atom (..))+import           Data.Map            (Map)+import qualified Data.Map as M       (fromList, (!))+import           Data.List           (groupBy, sortOn)+import           Data.Vector         (Vector)+import qualified Data.Vector as V    (toList)++groupChainByResidue :: Vector PDB.Atom -> [[PDB.Atom]]+groupChainByResidue = sortOn (sortOnResidue . head) . groupBy atomsFromSameResidue . V.toList+  where +    atomsFromSameResidue :: PDB.Atom -> PDB.Atom -> Bool+    atomsFromSameResidue atom1 atom2 = PDB.atomResSeq atom1 == PDB.atomResSeq atom2 && PDB.atomICode atom1 == PDB.atomICode atom2+    sortOnResidue :: PDB.Atom -> Int+    sortOnResidue PDB.Atom{..} = atomSerial * 100 + (insertionCodeSortingCorrections M.! atomICode)+    insertionCodeSortingCorrections :: Map Char Int+    insertionCodeSortingCorrections = M.fromList $ zip (' ':['A'..'Z']) [0..]

@@ -64,8 +64,13 @@  atomP :: Parser CoordLike atomP = let atom = Atom <$>-                    (string "ATOM " *>                                       -- (1 -  6) # we increased atomSerial field for one symbol-                    (read <$> count 6 notEndLineChar) <* space)              -- (7 - 11)  atomSerial+                    (+                      (string "ATOM " *>                                     -- (1 -  5)  ATOM -- we extended atomSerial length to the left for one symbol+                      (read <$> count 6 notEndLineChar) <* space)            -- (6 - 11)  atomSerial+                        <|>                                                  -- or+                      (string "HETATM" *>                                    -- (1 -  6)  HETATM+                      (read <$> count 5 notEndLineChar) <* space)            -- (7 - 11)  atomSerial+                    )                     <*> (T.pack <$> count 4 notEndLineChar)                  -- (13 - 16) atomName                     <*> notEndLineChar                                       -- (17)      atomAltLoc                     <*> (T.pack <$> count 3 notEndLineChar) <* space         -- (18 - 20) atomResName@@ -84,7 +89,7 @@  coordNotAtomP :: Parser CoordLike coordNotAtomP = do-    _ <- string "HETATM" <|> string "TER " <|> string "ANISOU" <|> string "CONECT"+    _ <- string "TER " <|> string "ANISOU" <|> string "CONECT"     skipWhile $ not . isEndOfLine     endOfLine     return CoordNotAtomLine

@@ -1,6 +1,7 @@ module Bio.PDB.Reader   ( fromTextPDB   , fromFilePDB+  , PDBWarnings(..)   ) where  import           Bio.PDB.Parser         (pdbP)@@ -61,13 +62,13 @@   then Right standardizedText   else Left "There are trash strings between model strings" -fromFilePDB :: MonadIO m => FilePath -> m (Either Text ([PDBWarnings], Either Text PDB))-fromFilePDB f = do-  content <- liftIO (TIO.readFile f)-  let preprocessed = preprocess content-  pure $ fmap (first pack . parseOnly pdbP) <$> preprocessed -fromTextPDB :: Text -> Either Text ([PDBWarnings], Either Text PDB)-fromTextPDB text = fmap (first pack . parseOnly pdbP) <$> preprocessed-  where-    preprocessed = preprocess text+fromFilePDB :: MonadIO m => FilePath -> m (Either Text ([PDBWarnings], PDB))+fromFilePDB = liftIO . fmap fromTextPDB . TIO.readFile++fromTextPDB :: Text -> Either Text ([PDBWarnings], PDB)+fromTextPDB text = do+  (warnings, preprocessedText) <- preprocess text+  pdb <- first pack $ parseOnly pdbP preprocessedText++  pure (warnings, pdb)

@@ -92,3 +92,6 @@                  , atomCharge     :: Text    -- ^ Charge on the atom.                  }   deriving (Show, Eq, Generic, NFData)++instance Ord Atom where+  a1 <= a2 = atomSerial a1 <= atomSerial a2

@@ -44,16 +44,20 @@  -- | Generic atom representation ---data Atom = Atom { atomId       :: GlobalID -- ^ global identifier-                 , atomName     :: Text     -- ^ IUPAC atom name-                 , atomElement  :: Text     -- ^ atom chemical element-                 , atomCoords   :: V3 Float -- ^ 3D coordinates of atom-                 , formalCharge :: Int      -- ^ Formal charge of atom-                 , bFactor      :: Float    -- ^ B-factor of atom-                 , occupancy    :: Float    -- ^ the amount of each conformation that is observed in the crystal+data Atom = Atom { atomId         :: GlobalID -- ^ global identifier, 0-based+                 , atomInputIndex :: Int      -- ^ atom index from input file+                 , atomName       :: Text     -- ^ IUPAC atom name+                 , atomElement    :: Text     -- ^ atom chemical element+                 , atomCoords     :: V3 Float -- ^ 3D coordinates of atom+                 , formalCharge   :: Int      -- ^ Formal charge of atom+                 , bFactor        :: Float    -- ^ B-factor of atom+                 , occupancy      :: Float    -- ^ the amount of each conformation that is observed in the crystal                  }   deriving (Show, Eq, Generic) +instance Ord Atom where+  a1 <= a2 = atomId a1 <= atomId a2+ instance NFData Atom  -- | Generic chemical bond@@ -76,11 +80,13 @@  -- | A set of atoms, organized to a residues ---data Residue = Residue { resName         :: Text                  -- ^ residue name-                       , resAtoms        :: Vector Atom           -- ^ a set of residue atoms-                       , resBonds        :: Vector (Bond LocalID) -- ^ a set of residue bonds with local identifiers (position in 'resAtoms')-                       , resSecondary    :: SecondaryStructure    -- ^ residue secondary structure-                       , resChemCompType :: Text                  -- ^ chemical component type+data Residue = Residue { resName          :: Text                  -- ^ residue name+                       , resNumber        :: Int                   -- ^ residue number+                       , resInsertionCode :: Char                  -- ^ residue insertion code+                       , resAtoms         :: Vector Atom           -- ^ a set of residue atoms+                       , resBonds         :: Vector (Bond LocalID) -- ^ a set of residue bonds with local identifiers (position in 'resAtoms')+                       , resSecondary     :: SecondaryStructure    -- ^ residue secondary structure+                       , resChemCompType  :: Text                  -- ^ chemical component type                        }   deriving (Show, Eq, Generic, NFData) 

@@ -5,86 +5,107 @@  module MAESpec where -import           Bio.MAE       (fromFile)+import           Bio.MAE       (modelsFromMaeFile) import           Bio.Structure (Atom (..), Bond (..), Chain (..), GlobalID (..),-                                LocalID (..), Model (..), Residue (..),-                                StructureModels (..))+                                LocalID (..), Model (..), Residue (..))++import           Control.Monad.IO.Class (MonadIO)+ import           Data.Set      (Set) import qualified Data.Set      as S (fromList) import           Data.Vector   (Vector)-import qualified Data.Vector   as V (fromList, toList, (!))+import qualified Data.Vector   as V (fromList, toList, (!), head)+import           Data.Either   (fromRight) import           Linear.V3     (V3 (..))+ import           Test.Hspec  maeSpec :: Spec-maeSpec = describe "Mae spec." $ do-    Model{..} <- runIO $ V.toList . modelsOf <$> fromFile "test/MAE/small.mae" >>= \[x] -> pure x+maeSpec = describe "Mae spec." $+    beforeAll (firstMaeModel "test/MAE/small.mae") $ do+        let firstChainResidues Model{..} = (chainResidues $ modelChains V.! 0)+        let secondChainResidues Model{..} = (chainResidues $ modelChains V.! 1) -    it "two chains" $ length modelChains `shouldBe` 2-    it "residues" $ do-        fmap resName (chainResidues $ modelChains V.! 0) `shouldBe` V.fromList ["ACE", "ASP", "ILE", "LYS"]-        fmap resName (chainResidues $ modelChains V.! 1) `shouldBe` V.fromList ["GLU", "LEU", "VAL", "ARG", "PRO", "GLY", "ALA", "LEU", "VAL"]+        it "two chains" $ \Model{..} -> length modelChains `shouldBe` 2+        +        it "residue numbers" $ \m@Model{..} -> do+            fmap resName (firstChainResidues m) `shouldBe` V.fromList ["ACE", "ASP", "ILE", "LYS"]+            fmap resName (secondChainResidues m) `shouldBe` V.fromList ["GLU", "LEU", "VAL", "ARG", "PRO", "GLY", "ALA", "LEU", "VAL"]+        +        it "residue names" $ \m@Model{..} -> do+            fmap resNumber (firstChainResidues m) `shouldBe` V.fromList [0, 1, 2, 3]+            fmap resNumber (secondChainResidues m) `shouldBe` V.fromList [10, 11, 12, 13, 13, 15, 16, 17, 18] -- 13 is doubled because the second 13 has 'A' insertion code -    it "atoms count" $ do-        sum (fmap (length . resAtoms) $ chainResidues $ modelChains V.! 0) `shouldBe` 62-        sum (fmap (length . resAtoms) $ chainResidues $ modelChains V.! 1) `shouldBe` 140+        it "residue insertion codes" $ \m@Model{..} -> do+            fmap resInsertionCode (firstChainResidues m) `shouldBe` V.fromList [' ', ' ', ' ', ' ']+            fmap resInsertionCode (secondChainResidues m) `shouldBe` V.fromList [' ', ' ', ' ', ' ', 'A', ' ', ' ', ' ', ' '] -    let allBonds = [ (2, 15, 1), (31, 44, 1), (15, 60, 1), (52, 27, 1), (57, 30, 1), (56, 8, 1), (38, 9, 1), (16, 36, 1), (31, 35, 1), (42, 22, 1), (38, 10, 1), (46, 40, 1), (36, 8, 1)-                   , (46, 7, 1), (35, 13, 1), (60, 31, 1), (40, 32, 1), (46, 45, 1), (29, 1, 1), (9, 53, 1), (22, 33, 1), (41, 18, 1), (45, 31, 1), (13, 53, 1), (40, 19, 1), (55, 20, 1)-                   , (16, 22, 1), (51, 36, 1), (46, 56, 1), (8, 17, 1), (56, 40, 1), (50, 6, 1), (57, 16, 1), (57, 34, 1), (48, 47, 1), (28, 48, 1), (57, 30, 1), (39, 40, 1), (11, 6, 1)-                   , (35, 55, 1), (47, 30, 1), (4, 35, 1), (60, 35, 1), (39, 56, 1), (44, 24, 1), (29, 55, 1), (29, 41, 1), (6, 55, 1), (52, 51, 1), (32, 21, 1), (55, 55, 1), (43, 55, 1)-                   , (30, 44, 1), (54, 47, 1), (13, 50, 1), (14, 56, 1), (44, 54, 1), (58, 21, 1), (27, 58, 1), (84, 122, 2), (122, 122, 2), (109, 79, 2), (103, 121, 2), (102, 139, 2)-                   , (93, 93, 2), (68, 103, 2), (110, 77, 2), (109, 63, 2), (86, 83, 2), (123, 137, 2), (110, 131, 2), (122, 85, 2), (75, 110, 2), (131, 138, 2), (88, 134, 2), (117, 81, 2)-                   , (107, 74, 2), (67, 61, 2), (97, 134, 2), (94, 131, 2), (65, 95, 2), (124, 100, 2), (120, 106, 2), (71, 111, 2), (95, 129, 2), (104, 116, 2), (61, 95, 2), (88, 132, 2)-                   , (80, 97, 2), (68, 121, 2), (138, 87, 2), (84, 134, 2), (86, 139, 2), (84, 76, 2), (83, 85, 2), (76, 61, 2), (116, 98, 2), (64, 82, 2), (106, 93, 2), (96, 102, 2), (98, 122, 2)-                   , (74, 82, 2), (123, 130, 2), (114, 127, 2), (102, 122, 2), (90, 126, 2), (118, 92, 2), (88, 71, 2), (120, 132, 2), (135, 71, 2), (125, 136, 2), (71, 67, 2), (128, 110, 2)-                   , (67, 95, 2), (83, 72, 2), (60, 70, 2), (99, 135, 2), (87, 87, 2), (63, 124, 2), (64, 110, 2), (127, 81, 2), (90, 106, 2), (93, 92, 2), (62, 91, 2), (67, 119, 2), (95, 70, 2)-                   , (134, 115, 2), (113, 115, 2), (93, 112, 2), (102, 106, 2), (136, 105, 2), (82, 126, 2), (109, 103, 2), (66, 119, 2), (84, 60, 2), (79, 67, 2), (63, 96, 2), (74, 134, 2), (112, 91, 2)-                   , (105, 68, 2), (125, 80, 2), (136, 75, 2), (62, 117, 2), (66, 124, 2), (105, 120, 2), (141, 63, 2), (63, 90, 2), (139, 60, 2), (77, 73, 2), (81, 118, 2), (81, 72, 2), (123, 126, 2)-                   , (128, 136, 2), (83, 91, 2), (86, 124, 2), (140, 108, 2), (126, 98, 2), (67, 92, 2), (131, 69, 2), (69, 94, 2), (137, 66, 2), (99, 87, 2), (134, 134, 2), (134, 72, 2), (91, 138, 2)-                   , (89, 97, 2), (109, 128, 2), (123, 100, 2), (61, 64, 2), (68, 114, 2), (119, 115, 2), (122, 106, 2), (131, 84, 2), (93, 94, 2), (124, 118, 2), (82, 119, 2), (107, 95, 2), (112, 64, 2)-                   , (71, 71, 2), (89, 62, 2), (89, 61, 2), (120, 80, 2), (85, 61, 2), (97, 95, 2), (94, 98, 2), (98, 139, 2), (115, 133, 2), (115, 64, 2), (136, 66, 2), (131, 62, 2), (83, 102, 2)-                   , (121, 60, 2), (117, 110, 2), (134, 138, 2), (119, 89, 2), (137, 128, 2), (112, 118, 2), (105, 116, 2), (95, 69, 2)-                   ]+        it "atoms count" $ \m@Model{..} -> do+            sum (fmap (length . resAtoms) (firstChainResidues m)) `shouldBe` 62+            sum (fmap (length . resAtoms) (secondChainResidues m)) `shouldBe` 140 -    it "global bonds" $ do-        modelBonds `shouldBe` toBond GlobalID <$> V.fromList allBonds+        let allBonds = [ (2, 15, 1), (31, 44, 1), (15, 60, 1), (52, 27, 1), (57, 30, 1), (56, 8, 1), (38, 9, 1), (16, 36, 1), (31, 35, 1), (42, 22, 1), (38, 10, 1), (46, 40, 1), (36, 8, 1)+                       , (46, 7, 1), (35, 13, 1), (60, 31, 1), (40, 32, 1), (46, 45, 1), (29, 1, 1), (9, 53, 1), (22, 33, 1), (41, 18, 1), (45, 31, 1), (13, 53, 1), (40, 19, 1), (55, 20, 1)+                       , (16, 22, 1), (51, 36, 1), (46, 56, 1), (8, 17, 1), (56, 40, 1), (50, 6, 1), (57, 16, 1), (57, 34, 1), (48, 47, 1), (28, 48, 1), (57, 30, 1), (39, 40, 1), (11, 6, 1)+                       , (35, 55, 1), (47, 30, 1), (4, 35, 1), (60, 35, 1), (39, 56, 1), (44, 24, 1), (29, 55, 1), (29, 41, 1), (6, 55, 1), (52, 51, 1), (32, 21, 1), (55, 55, 1), (43, 55, 1)+                       , (30, 44, 1), (54, 47, 1), (13, 50, 1), (14, 56, 1), (44, 54, 1), (58, 21, 1), (27, 58, 1), (84, 122, 2), (122, 122, 2), (109, 79, 2), (103, 121, 2), (102, 139, 2)+                       , (93, 93, 2), (68, 103, 2), (110, 77, 2), (109, 63, 2), (86, 83, 2), (123, 137, 2), (110, 131, 2), (122, 85, 2), (75, 110, 2), (131, 138, 2), (88, 134, 2), (117, 81, 2)+                       , (107, 74, 2), (67, 61, 2), (97, 134, 2), (94, 131, 2), (65, 95, 2), (124, 100, 2), (120, 106, 2), (71, 111, 2), (95, 129, 2), (104, 116, 2), (61, 95, 2), (88, 132, 2)+                       , (80, 97, 2), (68, 121, 2), (138, 87, 2), (84, 134, 2), (86, 139, 2), (84, 76, 2), (83, 85, 2), (76, 61, 2), (116, 98, 2), (64, 82, 2), (106, 93, 2), (96, 102, 2), (98, 122, 2)+                       , (74, 82, 2), (123, 130, 2), (114, 127, 2), (102, 122, 2), (90, 126, 2), (118, 92, 2), (88, 71, 2), (120, 132, 2), (135, 71, 2), (125, 136, 2), (71, 67, 2), (128, 110, 2)+                       , (67, 95, 2), (83, 72, 2), (60, 70, 2), (99, 135, 2), (87, 87, 2), (63, 124, 2), (64, 110, 2), (127, 81, 2), (90, 106, 2), (93, 92, 2), (62, 91, 2), (67, 119, 2), (95, 70, 2)+                       , (134, 115, 2), (113, 115, 2), (93, 112, 2), (102, 106, 2), (136, 105, 2), (82, 126, 2), (109, 103, 2), (66, 119, 2), (84, 60, 2), (79, 67, 2), (63, 96, 2), (74, 134, 2), (112, 91, 2)+                       , (105, 68, 2), (125, 80, 2), (136, 75, 2), (62, 117, 2), (66, 124, 2), (105, 120, 2), (141, 63, 2), (63, 90, 2), (139, 60, 2), (77, 73, 2), (81, 118, 2), (81, 72, 2), (123, 126, 2)+                       , (128, 136, 2), (83, 91, 2), (86, 124, 2), (140, 108, 2), (126, 98, 2), (67, 92, 2), (131, 69, 2), (69, 94, 2), (137, 66, 2), (99, 87, 2), (134, 134, 2), (134, 72, 2), (91, 138, 2)+                       , (89, 97, 2), (109, 128, 2), (123, 100, 2), (61, 64, 2), (68, 114, 2), (119, 115, 2), (122, 106, 2), (131, 84, 2), (93, 94, 2), (124, 118, 2), (82, 119, 2), (107, 95, 2), (112, 64, 2)+                       , (71, 71, 2), (89, 62, 2), (89, 61, 2), (120, 80, 2), (85, 61, 2), (97, 95, 2), (94, 98, 2), (98, 139, 2), (115, 133, 2), (115, 64, 2), (136, 66, 2), (131, 62, 2), (83, 102, 2)+                       , (121, 60, 2), (117, 110, 2), (134, 138, 2), (119, 89, 2), (137, 128, 2), (112, 118, 2), (105, 116, 2), (95, 69, 2)+                       ] -    let residue1 = (chainResidues $ modelChains V.! 0) V.! 3-    let residue2 = (chainResidues $ modelChains V.! 1) V.! 0+        it "global bonds" $ \Model{..} ->+            modelBonds `shouldBe` toBond GlobalID <$> V.fromList allBonds -    it "atoms in residues" $ do-        let atoms1 = resAtoms residue1-        let atoms2 = resAtoms residue2+        let residue1 m = firstChainResidues m V.! 3+        let residue2 m = secondChainResidues m V.! 0 -        fmap atomId atoms1 `shouldBe` GlobalID <$> V.fromList ([37 .. 58] <> [199 .. 201])-        fmap atomId atoms2 `shouldBe` GlobalID <$> V.fromList [59 .. 73]+        it "atoms in residues" $ \m@Model{..} -> do+            let atoms1 = resAtoms $ residue1 m+            let atoms2 = resAtoms $ residue2 m -        fmap atomName atoms1 `shouldBe` V.fromList [ "N", "CA", "C", "O", "CB", "CG", "CD", "CE", "NZ", "H"-                                                   , "HA", "HB3", "HB2", "HG3", "HG2", "HD3", "HD2", "HE3"-                                                   , "HE2", "HZ1", "HZ2", "HZ3", "2H", "CG2", "HA"-                                                   ]-        fmap atomName atoms2 `shouldBe` V.fromList [ "N", "CA", "C", "O", "CB", "CG", "CD", "OE1", "OE2", "H"-                                                   , "HA", "HB3", "HB2", "HG3", "HG2"-                                                   ]+            fmap atomId atoms1 `shouldBe` GlobalID <$> V.fromList ([37 .. 58] <> [199 .. 201])+            fmap atomId atoms2 `shouldBe` GlobalID <$> V.fromList [59 .. 73] -        fmap atomCoords atoms1 `shouldBe` V.fromList [ V3 (-3.444000) 7.750000 48.589000, V3 (-2.336000) 8.349000 47.850000, V3 (-1.127000) 7.408000 47.751000, V3 (-0.731000) 6.818000 48.755000, V3 (-1.974000) 9.706000 48.503000-                                                     , V3 (-0.665000) 10.349000 47.988000, V3 (-0.475000) 11.818000 48.390000, V3 (-1.315000) 12.784000 47.541000, V3 (-1.050000) 14.186000 47.903000, V3 (-3.565000) 8.072000 49.540000-                                                     , V3 (-2.692000) 8.558000 46.839000, V3 (-1.882000) 9.581000 49.583000, V3 (-2.810000) 10.388000 48.353000, V3 (-0.610000) 10.265000 46.901000, V3 0.184000 9.784000 48.375000-                                                     , V3 0.582000 12.073000 48.302000, V3 (-0.729000) 11.939000 49.444000, V3 (-2.380000) 12.585000 47.665000, V3 (-1.084000) 12.653000 46.483000, V3 (-1.609000) 14.795000 47.322000-                                                     , V3 (-1.289000) 14.335000 48.873000, V3 (-0.072000) 14.393000 47.762000, V3 (-9.600000) 7.518000 44.746000, V3 (-6.594000) 4.327000 47.986000, V3 (-5.134000) 6.683000 50.072000-                                                     ]-        fmap atomCoords atoms2 `shouldBe` V.fromList [ V3 (-8.382000) 11.633000 16.946000, V3 (-9.715000) 12.157000 17.191000, V3 (-9.590000) 13.665000 17.450000, V3 (-9.030000) 14.054000 18.475000, V3 (-10.323000) 11.388000 18.396000-                                                     , V3 (-11.833000) 11.576000 18.604000, V3 (-12.658000) 10.923000 17.495000, V3 (-12.525000) 9.689000 17.338000, V3 (-13.405000) 11.669000 16.827000, V3 (-7.615000) 12.187000 17.300000-                                                     , V3 (-10.322000) 11.996000 16.297000, V3 (-9.798000) 11.676000 19.305000, V3 (-10.128000) 10.318000 18.313000, V3 (-12.084000) 12.633000 18.684000, V3 (-12.123000) 11.111000 19.544000-                                                     ]+            fmap atomName atoms1 `shouldBe` V.fromList [ "N", "CA", "C", "O", "CB", "CG", "CD", "CE", "NZ", "H"+                                                       , "HA", "HB3", "HB2", "HG3", "HG2", "HD3", "HD2", "HE3"+                                                       , "HE2", "HZ1", "HZ2", "HZ3", "2H", "CG2", "HA"+                                                       ]+            fmap atomName atoms2 `shouldBe` V.fromList [ "N", "CA", "C", "O", "CB", "CG", "CD", "OE1", "OE2", "H"+                                                       , "HA", "HB3", "HB2", "HG3", "HG2"+                                                       ] -    it "bonds in residues" $ do-        toSet (resBonds residue1) `shouldBe` toSet (V.fromList $ toBond LocalID . (\(x, y, o) -> (x - 37, y - 37, o)) <$> filter (\(x, y, _) -> x - 1 `elem` [37 .. 58] && y - 1 `elem` [37 .. 58]) allBonds)-        toSet (resBonds residue2) `shouldBe` toSet (V.fromList $ toBond LocalID . (\(x, y, o) -> (x - 59, y - 59, o)) <$> filter (\(x, y, _) -> x - 1 `elem` [59 .. 73] && y - 1 `elem` [59 .. 73]) allBonds)-  where-    toBond :: (Int -> b) -> (Int, Int, Int) -> Bond b-    toBond f (x, y, o) = Bond (f $ x - 1) (f $ y - 1) o+            fmap atomCoords atoms1 `shouldBe` V.fromList [ V3 (-3.444000) 7.750000 48.589000, V3 (-2.336000) 8.349000 47.850000, V3 (-1.127000) 7.408000 47.751000, V3 (-0.731000) 6.818000 48.755000, V3 (-1.974000) 9.706000 48.503000+                                                         , V3 (-0.665000) 10.349000 47.988000, V3 (-0.475000) 11.818000 48.390000, V3 (-1.315000) 12.784000 47.541000, V3 (-1.050000) 14.186000 47.903000, V3 (-3.565000) 8.072000 49.540000+                                                         , V3 (-2.692000) 8.558000 46.839000, V3 (-1.882000) 9.581000 49.583000, V3 (-2.810000) 10.388000 48.353000, V3 (-0.610000) 10.265000 46.901000, V3 0.184000 9.784000 48.375000+                                                         , V3 0.582000 12.073000 48.302000, V3 (-0.729000) 11.939000 49.444000, V3 (-2.380000) 12.585000 47.665000, V3 (-1.084000) 12.653000 46.483000, V3 (-1.609000) 14.795000 47.322000+                                                         , V3 (-1.289000) 14.335000 48.873000, V3 (-0.072000) 14.393000 47.762000, V3 (-9.600000) 7.518000 44.746000, V3 (-6.594000) 4.327000 47.986000, V3 (-5.134000) 6.683000 50.072000+                                                         ]+            fmap atomCoords atoms2 `shouldBe` V.fromList [ V3 (-8.382000) 11.633000 16.946000, V3 (-9.715000) 12.157000 17.191000, V3 (-9.590000) 13.665000 17.450000, V3 (-9.030000) 14.054000 18.475000, V3 (-10.323000) 11.388000 18.396000+                                                         , V3 (-11.833000) 11.576000 18.604000, V3 (-12.658000) 10.923000 17.495000, V3 (-12.525000) 9.689000 17.338000, V3 (-13.405000) 11.669000 16.827000, V3 (-7.615000) 12.187000 17.300000+                                                         , V3 (-10.322000) 11.996000 16.297000, V3 (-9.798000) 11.676000 19.305000, V3 (-10.128000) 10.318000 18.313000, V3 (-12.084000) 12.633000 18.684000, V3 (-12.123000) 11.111000 19.544000+                                                         ] -    toSet :: Ord a => Vector a -> Set a-    toSet = S.fromList . V.toList+        it "bonds in residues" $ \m@Model{..} -> do+            toSet (resBonds $ residue1 m) `shouldBe` toSet (V.fromList $ toBond LocalID . (\(x, y, o) -> (x - 37, y - 37, o)) <$> filter (\(x, y, _) -> x - 1 `elem` [37 .. 58] && y - 1 `elem` [37 .. 58]) allBonds)+            toSet (resBonds $ residue2 m) `shouldBe` toSet (V.fromList $ toBond LocalID . (\(x, y, o) -> (x - 59, y - 59, o)) <$> filter (\(x, y, _) -> x - 1 `elem` [59 .. 73] && y - 1 `elem` [59 .. 73]) allBonds)+    where+        toBond :: (Int -> b) -> (Int, Int, Int) -> Bond b+        toBond f (x, y, o) = Bond (f $ x - 1) (f $ y - 1) o++        toSet :: Ord a => Vector a -> Set a+        toSet = S.fromList . V.toList++firstMaeModel :: (MonadIO m) => FilePath -> m Model+firstMaeModel filepath = do+  eitherMae <- modelsFromMaeFile filepath+  -- `evaluate . force` fails for some reason+  pure . V.head $ fromRight undefined eitherMae

@@ -2,8 +2,9 @@  import           Bio.MMTF import           Bio.MMTF.Decode.Codec-import           Data.Int              (Int8)-import           Data.Vector           ((!))+import           Data.Int               (Int8)+import           Data.Vector            ((!))+import           Data.ByteString.Lazy as BS (readFile) import           Test.Hspec  mmtfCodecSpec :: Spec@@ -26,11 +27,10 @@ mmtfParserSpec =   describe "MMTF parser" $   it "should parse 1FSD" $ do-    m <- fetch "1FSD"+    m <- BS.readFile "test/MMTF/1FSD.dms" >>= decode     (structureId . structure) m `shouldBe` "1FSD"     (numModels . structure) m `shouldBe` 41     (length . bFactorList . atom) m `shouldBe` 20664     ((! 0) . experimentalMethods . structure) m `shouldBe` "SOLUTION NMR"     ((! 0) . xCoordList . atom) m `shouldBe` (-12.847)     ((! 20663) . xCoordList . atom) m `shouldBe` 5.672-

@@ -0,0 +1,446 @@+{-# LANGUAGE OverloadedStrings #-}++module PDBParserSpec where++import           Bio.PDB.Reader  (fromTextPDB)+import           Bio.PDB.Type    (Atom (..), FieldType (..), PDB (..))+import qualified Data.Map.Strict (empty, fromList, singleton)+import           Data.Text       as T (Text, intercalate, length, lines, pack,+                                       replicate, take)+import qualified Data.Vector     as V (empty, fromList, singleton)+import           Test.Hspec++oneModelSpecP :: Spec+oneModelSpecP = describe "One model." $+        it "correctly parses pdb with only one model without strings \"MODEL\" & \"ENDMDL\"" $ do+        let mt  = fromTextPDB . lenghtenLines $ T.pack ( "HEADER header\n" +++                                         "TITLE STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME\n" +++                                         "COMPND compnd\n" +++                                         "SOURCE source\n" +++                                         "KEYWDS keywds\n" +++                                         "AUTHOR  Masha\n" +++                                         "REVDAT revdat\n" +++                                         "REMARK   1 REFERENCE 1\n" +++                                         "SEQRES seqres\n" +++                                         "CRYST1 cryst1\n" +++                                         "ORIGX1 origx1 n=1\n" +++                                         "SCALE2 sclaen n=2\n" +++                                         "ATOM   2032  OXT CYS A 214      -4.546 -29.673  26.796  1.0 143.51           O  \n" +++                                         "ATOM   2033  H   CYS A 214      -6.124 -27.225  26.558  1.00 15.00           H  \n" +++                                         "TER    2034      CYS A 214                                                      \n" +++                                         "ATOM   2035  N   GLU B   1      18.637 -61.583  66.852  1.0 118.48           N  \n" +++                                         "TER   12534      ARG D 474                                                      \n" +++                                         "HETATM12535  C1  NAG B 475       5.791 -20.194  -7.051  1.00 34.66           C  \n" +++                                         "HETATM12538  C4  NAG B 475       6.943 -19.507  -9.597  1.00 25.87           C  \n" +++                                         "ATOM   2036  CA  GLU B   1      19.722 -62.606  66.868  1.00 19.77           C  \n" +++                                         "ATOM   2036  CA  GLU B   1A     19.722 -62.606  66.868  1.00 19.77           C  \n" +++                                         "CONECT conect\n" +++                                         "CONECT conect conect\n" +++                                         "MASTER 1 2 3 4 5 6 7 8\n" +++                                         "END"+                                        )+        mt `shouldBe` Right ([], oneModelPDB)++oneModelPDB :: PDB+oneModelPDB =  PDB { title = "STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME"+                   , remarks = Data.Map.Strict.singleton (Just 1) (V.singleton "REFERENCE 1")+                   , models = V.singleton $ V.fromList [V.fromList [+                                                                Atom {atomSerial = 2032, atomName = " OXT", atomAltLoc = ' ', atomResName = "CYS", atomChainID = 'A', atomResSeq = 214, atomICode = ' ',+                                                                        atomX = -4.546, atomY = -29.673, atomZ = 26.796, atomOccupancy = 1.0, atomTempFactor = 143.51, atomElement = " O", atomCharge = "  "},+                                                                Atom {atomSerial = 2033, atomName = " H  ", atomAltLoc = ' ', atomResName = "CYS", atomChainID = 'A', atomResSeq = 214, atomICode = ' ', +                                                                        atomX = -6.124, atomY = -27.225, atomZ = 26.558, atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}],+                                                        V.fromList [+                                                                Atom {atomSerial = 2035, atomName = " N  ", atomAltLoc = ' ', atomResName = "GLU", atomChainID = 'B', atomResSeq = 1, atomICode = ' ',+                                                                        atomX = 18.637, atomY = -61.583, atomZ = 66.852, atomOccupancy = 1.0, atomTempFactor = 118.48, atomElement = " N", atomCharge = "  "},+                                                                Atom {atomSerial = 12535, atomName = " C1 ", atomAltLoc = ' ', atomResName = "NAG", atomChainID = 'B', atomResSeq = 475, atomICode = ' ',+                                                                        atomX = 5.791, atomY = -20.194, atomZ = -7.051, atomOccupancy = 1.0, atomTempFactor = 34.66, atomElement = " C", atomCharge = "  "},+                                                                Atom {atomSerial = 12538, atomName = " C4 ", atomAltLoc = ' ', atomResName = "NAG", atomChainID = 'B', atomResSeq = 475, atomICode = ' ', +                                                                        atomX = 6.943, atomY = -19.507, atomZ = -9.597, atomOccupancy = 1.0, atomTempFactor = 25.87, atomElement = " C", atomCharge = "  "},+                                                                Atom {atomSerial = 2036, atomName = " CA ", atomAltLoc = ' ', atomResName = "GLU", atomChainID = 'B', atomResSeq = 1, atomICode = ' ',+                                                                        atomX = 19.722, atomY = -62.606, atomZ = 66.868, atomOccupancy = 1.0, atomTempFactor = 19.77, atomElement = " C", atomCharge = "  "},+                                                                Atom {atomSerial = 2036, atomName = " CA ", atomAltLoc = ' ', atomResName = "GLU", atomChainID = 'B', atomResSeq = 1, atomICode = 'A',+                                                                        atomX = 19.722, atomY = -62.606, atomZ = 66.868, atomOccupancy = 1.0, atomTempFactor = 19.77, atomElement = " C", atomCharge = "  "}]+                                                        ]+                  , otherFields = Data.Map.Strict.fromList [(HEADER, V.fromList [" header"]),(COMPND,V.fromList[" compnd"]),(SOURCE,V.fromList[" source"]),(KEYWDS,V.fromList[" keywds"]),(AUTHOR,V.fromList["  Masha"]),(REVDAT,V.fromList[" revdat"]),(SEQRES,V.fromList [" seqres"]),+                                                            (CRYST1,V.fromList [" cryst1"]),(ORIGX1,V.fromList[" origx1 n=1"]),(SCALE2,V.fromList [" sclaen n=2"]),(MASTER,V.fromList[" 1 2 3 4 5 6 7 8"])]+                   }++manyModelsSpecP :: Spec+manyModelsSpecP = describe "Some models." $+        it "correctly parses pdb with many models - they have strings \"MODEL\" & \"ENDMDL\" and text has disordered strings" $ do+        let mt = fromTextPDB . lenghtenLines $ T.pack ( "HEADER header\n" +++                                        "TITLE STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME\n" +++                                        "COMPND compnd\n" +++                                        "SOURCE source\n" +++                                        "REMARK 1   REFERENCE 1\n" +++                                        "REMARK 2   Reference 2_1/1\n" +++                                        "CRYST1 cryst1\n" +++                                        "ORIGX1 origx1 n=1\n" +++                                        "SCALE2 sclaen n=2\n" +++                                        "MODEL  4  \n" +++                                        "ATOM   2033  H   CYS A 214      -6.124 -27.225  26.558  1.00 15.00           H  \n" +++                                        "TER    2034      CYS A 214                                                      \n" +++                                        "ANISOU anisou\n" +++                                        "ENDMDL \n" +++                                        "MODEL  5 \n" +++                                        "ATOM   2035  N   GLU B   1      18.637-691.583  66.852  1.0 118.48           N  \n" +++                                        "ATOM  12531 HH12 ARG D 474      45.558 -39.551 -49.936  1.00 15.00           H  \n" +++                                        "ENDMDL \n" +++                                        "MODEL  6 \n" +++                                        "ATOM  12532 HH21 ARG D 474      47.457 -38.007 -47.445  1.00 15.00           H  \n" +++                                        "TER   12534      ARG D 474                                                      \n" +++                                        "ATOM  12533 HH22 ARG D 474      47.405 -39.268 -48.629  1.00 15.00           H  \n" +++                                        "HETATM12535  C1  NAG B 475       5.791 -20.194  -7.051  1.00 34.66           C  \n" +++                                        "HETATM12538  C4  NAG B 475       6.943 -19.507  -9.597  1.00 25.87           C  \n" +++                                        "CONECT conect conect\n" +++                                        "CONECT conect\n" +++                                        "ENDMDL\n" +++                                        "SEQRES seqres\n" +++                                        "KEYWDS keywds\n" +++                                        "EXPDTA expdta\n" +++                                        "AUTHOR  Masha\n" +++                                        "REVDAT revdat\n" +++                                        "MASTER 1 2 3 4 5 6 7 8\n"+                                      )+        mt `shouldBe` Right ([], manyModelsPDB)++manyModelsPDB :: PDB+manyModelsPDB = PDB { title = "STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME"+                    , remarks = Data.Map.Strict.fromList [(Just 1, V.singleton "REFERENCE 1"), (Just 2, V.singleton "Reference 2_1/1")]+                    , models = V.fromList [V.fromList  [V.fromList [+                                                                Atom {atomSerial = 2033, atomName = " H  ", atomAltLoc = ' ', atomResName = "CYS", atomChainID = 'A', atomResSeq = 214, atomICode = ' ', +                                                                        atomX = -6.124, atomY = -27.225, atomZ = 26.558, atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}]+                                                        ],+                                           V.fromList  [V.fromList [+                                                                Atom {atomSerial = 2035, atomName = " N  ", atomAltLoc = ' ', atomResName = "GLU", atomChainID = 'B', atomResSeq = 1, atomICode = ' ', +                                                                        atomX = 18.637, atomY = -691.583, atomZ = 66.852, atomOccupancy = 1.0, atomTempFactor = 118.48, atomElement = " N", atomCharge = "  "}],+                                                        V.fromList [+                                                                Atom {atomSerial = 12531, atomName = "HH12", atomAltLoc = ' ', atomResName = "ARG", atomChainID = 'D', atomResSeq = 474, atomICode = ' ', +                                                                        atomX = 45.558, atomY = -39.551, atomZ = -49.936, atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}]+                                                        ],+                                           V.fromList  [V.fromList [+                                                                Atom {atomSerial = 12532, atomName = "HH21", atomAltLoc = ' ', atomResName = "ARG", atomChainID = 'D', atomResSeq = 474, atomICode = ' ', +                                                                        atomX = 47.457, atomY = -38.007, atomZ = -47.445, atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "},+                                                                Atom {atomSerial = 12533, atomName = "HH22", atomAltLoc = ' ', atomResName = "ARG", atomChainID = 'D', atomResSeq = 474, atomICode = ' ', +                                                                        atomX = 47.405, atomY = -39.268, atomZ = -48.629, atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}],+                                                        V.fromList [+                                                                Atom {atomSerial = 12535, atomName = " C1 ", atomAltLoc = ' ', atomResName = "NAG", atomChainID = 'B', atomResSeq = 475, atomICode = ' ',+                                                                        atomX = 5.791, atomY = -20.194, atomZ = -7.051, atomOccupancy = 1.0, atomTempFactor = 34.66, atomElement = " C", atomCharge = "  "},+                                                                Atom {atomSerial = 12538, atomName = " C4 ", atomAltLoc = ' ', atomResName = "NAG", atomChainID = 'B', atomResSeq = 475, atomICode = ' ', +                                                                        atomX = 6.943, atomY = -19.507, atomZ = -9.597, atomOccupancy = 1.0, atomTempFactor = 25.87, atomElement = " C", atomCharge = "  "}]+                                                        ]+                                           ]+                    , otherFields = Data.Map.Strict.fromList [(HEADER, V.fromList [" header"]),(COMPND,V.fromList[" compnd"]),(SOURCE,V.fromList[" source"]),(KEYWDS,V.fromList[" keywds"]),(EXPDTA,V.fromList[" expdta"]),(AUTHOR,V.fromList["  Masha"]),(REVDAT,V.fromList[" revdat"]),(SEQRES,V.fromList [" seqres"]),+                                                              (CRYST1,V.fromList [" cryst1"]),(ORIGX1,V.fromList[" origx1 n=1"]),(SCALE2,V.fromList [" sclaen n=2"]),(MASTER,V.fromList[" 1 2 3 4 5 6 7 8"])]+                    }++noModelsSpecP :: Spec+noModelsSpecP = describe "No models." $+        it "correctly parses pdb without models (no ATOM, TER, HETATM strings)" $ do+        let mt = fromTextPDB . lenghtenLines $ T.pack ( "HEADER header\n" +++                                        "TITLE STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME\n" +++                                        "COMPND compnd\n" +++                                        "SOURCE source\n" +++                                        "KEYWDS keywds\n" +++                                        "EXPDTA expdta\n" +++                                        "AUTHOR  Masha\n" +++                                        "REVDAT revdat\n" +++                                        "REMARK 1   REFERENCE 1\n" +++                                        "SEQRES seqres\n" +++                                        "CRYST1 cryst1\n" +++                                        "ORIGX1 origx1 n=1\n" +++                                        "SCALE2 sclaen n=2\n" +++                                        "MASTER 1 2 3 4 5 6 7 8\n"+                                      )+        mt `shouldBe` Right ([], noModelsPDB)++noModelsPDB :: PDB+noModelsPDB = PDB { title = "STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME"+                  , remarks = Data.Map.Strict.singleton (Just 1) (V.singleton "REFERENCE 1")+                  , models = V.empty+                  , otherFields = Data.Map.Strict.fromList [(HEADER, V.fromList [" header"]),(COMPND,V.fromList[" compnd"]),(SOURCE,V.fromList[" source"]),(KEYWDS,V.fromList[" keywds"]),(EXPDTA,V.fromList[" expdta"]),(AUTHOR,V.fromList["  Masha"]),(REVDAT,V.fromList[" revdat"]),(SEQRES,V.fromList [" seqres"]),+                                                            (CRYST1,V.fromList [" cryst1"]),(ORIGX1,V.fromList[" origx1 n=1"]),(SCALE2,V.fromList [" sclaen n=2"]),(MASTER,V.fromList[" 1 2 3 4 5 6 7 8"])]+                  }++allFieldsModelSpecP :: Spec+allFieldsModelSpecP = describe "PDB with all strings." $+        it "correctly parses pdb with all types of string" $ do+        let mt = fromTextPDB . lenghtenLines $ T.pack ( "HEADER header\n" +++                                        "OBSLTE obslte\n" +++                                        "TITLE STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME\n" +++                                        "SPLIT split\n" +++                                        "CAVEAT caveat\n" +++                                        "COMPND compnd\n" +++                                        "SOURCE source\n" +++                                        "KEYWDS keywds\n" +++                                        "EXPDTA expdta1\n" +++                                        "EXPDTA expdta2\n" +++                                        "NUMMDL nummdl\n" +++                                        "MDLTYP mdltyp mdltyp\n" +++                                        "AUTHOR  Masha\n" +++                                        "REVDAT revdat\n" +++                                        "SPRSDE sprsde\n" +++                                        "JRNL   jrnl\n" +++                                        "REMARK 1   REFERENCE 1\n" +++                                        "DBREF dbref\n" +++                                        "DBREF1 dbref1\n" +++                                        "DBREF2 dbref2_1\n" +++                                        "DBREF2 dbref2_2\n" +++                                        "SEQADV seqadv\n" +++                                        "SEQRES seqres\n" +++                                        "MODRES modres\n" +++                                        "HET    het\n" +++                                        "HETNAM hetnam\n" +++                                        "HETSYN hetsyn\n" +++                                        "FORMUL of love\n" +++                                        "HELIX helix\n" +++                                        "SHEET sheet\n" +++                                        "SSBOND ssbond\n" +++                                        "LINK   link\n" +++                                        "CISPEP cispep\n" +++                                        "SITE   site\n" +++                                        "CRYST1 cryst1\n" +++                                        "ORIGX1 origx1 n=1\n" +++                                        "SCALE2 sclaen n=2\n" +++                                        "MTRIX3 matrixn n=3\n" +++                                        "MODEL 1\n" +++                                        "ATOM  12532 HH21 ARG D 474      47.457 -38.007 -47.445  1.00 15.00           H  \n" +++                                        "TER   12534      ARG D 474                                                      \n" +++                                        "HETATM12535  C1  NAG B 475       5.791 -20.194  -7.051  1.00 34.66           C  \n" +++                                        "CONECT conect\n" +++                                        "CONECT conect conect\n" +++                                        "ENDMDL\n" +++                                        "MASTER 1 2 3 4 5 6 7 8\n" +++                                        "END"+                                      )+        mt `shouldBe` Right ([], pdbWithAllFields)++pdbWithAllFields :: PDB+pdbWithAllFields = PDB { title = "STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME"+                       , models = V.fromList   [V.fromList [+                                                        V.singleton Atom {atomSerial = 12532, atomName = "HH21", atomAltLoc = ' ', atomResName = "ARG", atomChainID = 'D', atomResSeq = 474, atomICode = ' ', +                                                                                atomX = 47.457, atomY = -38.007, atomZ = -47.445, atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "},+                                                        V.singleton Atom {atomSerial = 12535, atomName = " C1 ", atomAltLoc = ' ', atomResName = "NAG", atomChainID = 'B', atomResSeq = 475, atomICode = ' ',+                                                                                atomX = 5.791, atomY = -20.194, atomZ = -7.051, atomOccupancy = 1.0, atomTempFactor = 34.66, atomElement = " C", atomCharge = "  "}]+                                                ]+                       , remarks = Data.Map.Strict.singleton (Just 1) (V.singleton "REFERENCE 1")+                       , otherFields = Data.Map.Strict.fromList [(HEADER, V.fromList [" header"]),(OBSLTE, V.fromList [" obslte"]), (SPLIT, V.fromList ["split"]), (CAVEAT, V.fromList [" caveat"]),+                                                                 (COMPND, V.fromList [" compnd"]),(SOURCE, V.fromList [" source"]),(KEYWDS, V.fromList [" keywds"]),+                                                                 (EXPDTA, V.fromList [" expdta1"," expdta2"]),(NUMMDL, V.fromList [" nummdl"]),(MDLTYP, V.fromList [" mdltyp mdltyp"]),+                                                                 (AUTHOR, V.fromList ["  Masha"]),(REVDAT, V.fromList [" revdat"]),(SPRSDE, V.fromList [" sprsde"]),+                                                                 (JRNL, V.fromList [" jrnl"]),(DBREF, V.fromList ["dbref"]), (DBREF1, V.fromList [" dbref1"]), (DBREF2, V.fromList [" dbref2_1"," dbref2_2"]),+                                                                 (SEQADV, V.fromList [" seqadv"]),(SEQRES, V.fromList [" seqres"]),(MODRES, V.fromList [" modres"]),+                                                                 (HET, V.fromList [" het"]), (HETNAM, V.fromList [" hetnam"]), (HETSYN, V.fromList [" hetsyn"]),(FORMUL, V.fromList [" of love"]),(HELIX, V.fromList ["helix"]),+                                                                 (SHEET, V.fromList ["sheet"]),(SSBOND, V.fromList [" ssbond"]),(LINK, V.fromList [" link"]),+                                                                 (CISPEP, V.fromList [" cispep"]),(SITE,V.fromList [" site"]),(CRYST1, V.fromList [" cryst1"]),(MTRIX3, V.fromList [" matrixn n=3"]),+                                                                 (ORIGX1, V.fromList [" origx1 n=1"]),(SCALE2, V.fromList [" sclaen n=2"]),(MASTER, V.fromList [" 1 2 3 4 5 6 7 8"])]+                       }+++emptySpecP :: Spec+emptySpecP = describe "empty PDB." $+        it "correctly parses empty pdb" $ do+        let mt = fromTextPDB ""+        mt `shouldBe` Right ([], emptyPdb)++emptyPdb :: PDB+emptyPdb = PDB { title = ""+               , models = V.empty+               , remarks = Data.Map.Strict.empty+               , otherFields = Data.Map.Strict.empty+               }++trashBetweenModelsSpecP :: Spec+trashBetweenModelsSpecP = describe "PDB has trash." $+        it "correctly parses pdb with trash string between models and other field strings" $ do+        let mt = fromTextPDB  . lenghtenLines $ T.pack ( "trash strings 1\n" +++                                        "HEADER header\n" +++                                        "TITLE STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME\n" +++                                        "REMARK 1   REFERENCE 1\n" +++                                        "trash strings 2\n" +++                                        "REMARK 2   Reference 2_1/1\n" +++                                        "MODEL  4  \n" +++                                        "ATOM   2033  H   CYS A 214      -6.124 -27.225  26.558  1.00 15.00           H  \n" +++                                        "TER    2034      CYS A 214                                                      \n" +++                                        "ANISOU anisou\n" +++                                        "ENDMDL \n" +++                                        "CRYST1 cryst1\n" +++                                        "trash strings 3\n" +++                                        "SEQRES seqres\n" +++                                        "MODEL  5 \n" +++                                        "ATOM   2035  N   GLU B   1      18.637-691.583  66.852  1.0 118.48           N  \n" +++                                        "ATOM  12531 HH12 ARG D 474      45.558 -39.551 -49.936  1.00 15.00           H  \n" +++                                        "ENDMDL \n" +++                                        "trash strings 4\n" +++                                        "MODEL  6 \n" +++                                        "ATOM  12532 HH21 ARG D 474      47.457 -38.007 -47.445  1.00 15.00           H  \n" +++                                        "TER   12534      ARG D 474                                                      \n" +++                                        "ATOM  12533 HH22 ARG D 474      47.405 -39.268 -48.629  1.00 15.00           H  \n" +++                                        "HETATM12535  C1  NAG B 475       5.791 -20.194  -7.051  1.00 34.66           C  \n" +++                                        "HETATM12538  C4  NAG B 475       6.943 -19.507  -9.597  1.00 25.87           C  \n" +++                                        "CONECT conect conect\n" +++                                        "ENDMDL\n" +++                                        "MASTER 1 2 3 4 5 6 7 8\n" +++                                        "trash strings 5\n"+                                       )+        mt `shouldBe` Left "There are trash strings between model strings"++pdbWithoutTrash :: PDB+pdbWithoutTrash = PDB { title = "STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME"+                      , remarks = Data.Map.Strict.fromList [(Just 1, V.singleton "REFERENCE 1"), (Just 2, V.singleton "Reference 2_1/1")]+                      , models = V.fromList [V.fromList [V.fromList [Atom {atomSerial = 2033, atomName = " H  ", atomAltLoc = ' ', atomResName = "CYS",+                                                                                                    atomChainID = 'A', atomResSeq = 214, atomICode = ' ', atomX = -6.124, atomY = -27.225,+                                                                                                    atomZ = 26.558, atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}]+                                                                            ],+                                                      V.fromList [V.fromList [Atom {atomSerial = 2035, atomName = " N  ", atomAltLoc = ' ', atomResName = "GLU",+                                                                                                   atomChainID = 'B', atomResSeq = 1, atomICode = ' ', atomX = 18.637, atomY = -691.583,+                                                                                                   atomZ = 66.852, atomOccupancy = 1.0, atomTempFactor = 118.48, atomElement = " N", atomCharge = "  "}],+                                                                            V.fromList [Atom {atomSerial = 12531, atomName = "HH12", atomAltLoc = ' ', atomResName = "ARG",+                                                                                                   atomChainID = 'D', atomResSeq = 474, atomICode = ' ', atomX = 45.558, atomY = -39.551, atomZ = -49.936,+                                                                                                   atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}]+                                                                           ],+                                                      V.fromList [V.fromList [Atom {atomSerial = 12532, atomName = "HH21", atomAltLoc = ' ', atomResName = "ARG",+                                                                                                   atomChainID = 'D', atomResSeq = 474, atomICode = ' ', atomX = 47.457, atomY = -38.007, atomZ = -47.445,+                                                                                                   atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "},+                                                                                                  Atom {atomSerial = 12533, atomName = "HH22", atomAltLoc = ' ', atomResName = "ARG",+                                                                                                   atomChainID = 'D', atomResSeq = 474, atomICode = ' ', atomX = 47.405, atomY = -39.268, atomZ = -48.629,+                                                                                                   atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}]+                                                                            ]+                                                        ]+                      , otherFields = Data.Map.Strict.fromList [(HEADER, V.fromList [" header"]),(SEQRES,V.fromList [" seqres"]),+                                                                (CRYST1,V.fromList [" cryst1"]),(MASTER,V.fromList[" 1 2 3 4 5 6 7 8"])]+                      }++onlyOneModelSpecP :: Spec+onlyOneModelSpecP = describe "Only One model." $+        it "correctly parses pdb with only one model without other field/title/trash/remarks strings" $ do+        let mt = fromTextPDB . lenghtenLines $ T.pack ( "ATOM   2032  OXT CYS A 214      -4.546 -29.673  26.796  1.0 143.51           O  \n" +++                                         "ATOM   2033  H   CYS A 214      -6.124 -27.225  26.558  1.00 15.00           H  \n" +++                                         "TER    2034      CYS A 214                                                      \n" +++                                         "ATOM   2035  N   GLU B   1      18.637 -61.583  66.852  1.0 118.48           N  \n" +++                                         "TER   12534      ARG D 474                                                      \n" +++                                         "HETATM12535  C1  NAG B 475       5.791 -20.194  -7.051  1.00 34.66           C  \n" +++                                         "HETATM12538  C4  NAG B 475       6.943 -19.507  -9.597  1.00 25.87           C  \n" +++                                         "ATOM   2036  CA  GLU B   1      19.722 -62.606  66.868  1.00 19.77           C  \n"+                                        )+        mt `shouldBe` Right ([], onlyOneModelPDB)++onlyOneModelPDB :: PDB+onlyOneModelPDB = PDB { title = ""+                      , remarks = Data.Map.Strict.empty+                      , models = V.singleton $ V.fromList [ +                                V.fromList [+                                        Atom {atomSerial = 2032, atomName = " OXT", atomAltLoc = ' ', atomResName = "CYS", atomChainID = 'A', atomResSeq = 214, atomICode = ' ',+                                                atomX = -4.546, atomY = -29.673, atomZ = 26.796, atomOccupancy = 1.0, atomTempFactor = 143.51, atomElement = " O", atomCharge = "  "},+                                        Atom {atomSerial = 2033, atomName = " H  ", atomAltLoc = ' ', atomResName = "CYS", atomChainID = 'A', atomResSeq = 214, atomICode = ' ', +                                                atomX = -6.124, atomY = -27.225, atomZ = 26.558, atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}],+                                V.fromList [+                                        Atom {atomSerial = 2035, atomName = " N  ", atomAltLoc = ' ', atomResName = "GLU", atomChainID = 'B', atomResSeq = 1, atomICode = ' ',+                                                atomX = 18.637, atomY = -61.583, atomZ = 66.852, atomOccupancy = 1.0, atomTempFactor = 118.48, atomElement = " N", atomCharge = "  "},+                                        Atom {atomSerial = 12535, atomName = " C1 ", atomAltLoc = ' ', atomResName = "NAG", atomChainID = 'B', atomResSeq = 475, atomICode = ' ',+                                                atomX = 5.791, atomY = -20.194, atomZ = -7.051, atomOccupancy = 1.0, atomTempFactor = 34.66, atomElement = " C", atomCharge = "  "},+                                        Atom {atomSerial = 12538, atomName = " C4 ", atomAltLoc = ' ', atomResName = "NAG", atomChainID = 'B', atomResSeq = 475, atomICode = ' ', +                                                atomX = 6.943, atomY = -19.507, atomZ = -9.597, atomOccupancy = 1.0, atomTempFactor = 25.87, atomElement = " C", atomCharge = "  "},+                                        Atom {atomSerial = 2036, atomName = " CA ", atomAltLoc = ' ', atomResName = "GLU", atomChainID = 'B', atomResSeq = 1, atomICode = ' ',+                                                atomX = 19.722, atomY = -62.606, atomZ = 66.868, atomOccupancy = 1.0, atomTempFactor = 19.77, atomElement = " C", atomCharge = "  "}]+                                                ]+                      , otherFields = Data.Map.Strict.empty+                      }++repeatedStringsSpecP :: Spec+repeatedStringsSpecP = describe "PDB with repeated other field strings." $+        it "correctly parses pdb with repeated other field strings (SOURCE)" $ do+        let mt = fromTextPDB . lenghtenLines $ T.pack ( "HEADER header\n" +++                                        "TITLE STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME\n" +++                                        "COMPND compnd\n" +++                                        "EXPDTA expdta\n" +++                                        "AUTHOR  Masha\n" +++                                        "SOURCE source1\n" +++                                        "KEYWDS keywds\n" +++                                        "REVDAT revdat\n" +++                                        "REMARK 1   REFERENCE 1\n" +++                                        "SEQRES seqres\n" +++                                        "CRYST1 cryst1\n" +++                                        "ORIGX1 origx1 n=1\n" +++                                        "SOURCE source2\n" +++                                        "SCALE2 sclaen n=2\n" +++                                        "MASTER 1 2 3 4 5 6 7 8\n" +++                                        "SOURCE source3\n" +++                                        "END"+                                      )+        mt `shouldBe` Right ([], repeatedStringsPDB)++repeatedStringsPDB :: PDB+repeatedStringsPDB = PDB { title = "STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME"+                         , remarks = Data.Map.Strict.singleton (Just 1) (V.singleton "REFERENCE 1")+                         , models = V.empty+                         , otherFields = Data.Map.Strict.fromList [(HEADER, V.fromList [" header"]),(COMPND,V.fromList[" compnd"]),(SOURCE,V.fromList[" source1", " source2", " source3"]),(KEYWDS,V.fromList[" keywds"]),(EXPDTA,V.fromList[" expdta"]),(AUTHOR,V.fromList["  Masha"]),(REVDAT,V.fromList[" revdat"]),(SEQRES,V.fromList [" seqres"]),+                                                                   (CRYST1,V.fromList [" cryst1"]),(ORIGX1,V.fromList[" origx1 n=1"]),(SCALE2,V.fromList [" sclaen n=2"]),(MASTER,V.fromList[" 1 2 3 4 5 6 7 8"])]+                         }++emptyRemarkSpecP :: Spec+emptyRemarkSpecP = describe "PDB with repeated remark strings without code." $+        it "correctly parses pdb with repeated remark strings without code" $ do+        let mt = fromTextPDB . lenghtenLines $ T.pack ( "REMARK 111 remark111_1/2\n" +++                                        "COMPND compnd\n" +++                                        "SOURCE source1\n" +++                                        "KEYWDS keywds\n" +++                                        "REMARK 111 remark111_2/2\n" +++                                        "REVDAT revdat\n" +++                                        "REMARK 2   remark2_1/1\n" +++                                        "SOURCE source2\n" +++                                        "REMARK     empty remark 1/2\n" +++                                        "END    \n" +++                                        "REMARK     empty remark 2/2\n" +++                                        "SCALE2 sclaen n=2\n"+                                      )+        mt `shouldBe` Right ([], pdbWithEmptyRemarks)++pdbWithEmptyRemarks :: PDB+pdbWithEmptyRemarks = PDB { title = ""+                         , remarks = Data.Map.Strict.fromList [(Just 2, V.fromList ["remark2_1/1"]), (Just 111, V.fromList ["remark111_1/2", "remark111_2/2"]), (Nothing, V.fromList ["empty remark 1/2", "empty remark 2/2"])]+                         , models = V.empty+                         , otherFields = Data.Map.Strict.fromList [(COMPND,V.fromList[" compnd"]),(SOURCE,V.fromList[" source1", " source2"]),+                                                                    (KEYWDS,V.fromList[" keywds"]), (REVDAT,V.fromList[" revdat"]), (SCALE2,V.fromList [" sclaen n=2"])]+                         }++emptyModelSpecP :: Spec+emptyModelSpecP = describe "PDB with one empty model." $+        it "correctly parses pdb with one model without strings inside" $ do+        let mt = fromTextPDB . lenghtenLines $ T.pack ( "trash strings 1\n" +++                                        "HEADER header\n" +++                                        "TITLE STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME\n" +++                                        "REMARK 1   REFERENCE 1\n" +++                                        "trash strings 2\n" +++                                        "REMARK 2   Reference 2_1/1\n" +++                                        "MODEL  4  \n" +++                                        "ENDMDL\n" +++                                        "MASTER 1 2 3 4 5 6 7 8\n" +++                                        "trash strings 5\n"+                                      )+        mt `shouldBe` Right ([], pdbEmptyModel)+++pdbEmptyModel :: PDB+pdbEmptyModel = PDB { title = "STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME"+                    , remarks = Data.Map.Strict.fromList [(Just 1, V.singleton "REFERENCE 1"), (Just 2, V.singleton "Reference 2_1/1")]+                    , models = V.empty+                    , otherFields = Data.Map.Strict.fromList [(HEADER, V.fromList [" header"]),+                                                              (MASTER,V.fromList[" 1 2 3 4 5 6 7 8"])]+                    }++lenghtenLines :: Text -> Text+lenghtenLines text = longLinedText+   where+       textLines = T.lines text+       longTextLines = changeLine <$> textLines+       desiredLength = 80  -- cause it is max length in standart pdb file+       longLinedText = T.intercalate "\n" longTextLines++       changeLine :: Text -> Text+       changeLine line | T.length line > desiredLength = T.take desiredLength line+                       | T.length line < desiredLength = line <> T.replicate spacesCount " "+                       | otherwise = line+            where+                    spacesCount = desiredLength - T.length line+

@@ -1,428 +1,149 @@ {-# LANGUAGE OverloadedStrings #-}+{-# LANGUAGE RecordWildCards   #-}  module PDBSpec where -import           Bio.PDB.Reader  (fromTextPDB)-import           Bio.PDB.Type    (Atom (..), FieldType (..), PDB (..))-import qualified Data.Map.Strict (empty, fromList, singleton)-import           Data.Text       as T (Text, intercalate, length, lines, pack,-                                       replicate, take)-import qualified Data.Vector     as V (empty, fromList, singleton)-import           Test.Hspec+import           Bio.PDB                (modelsFromPDBFile)+import           Bio.MAE                (modelsFromMaeFile)+import           Bio.Structure          ( Model(..), Chain(..), Residue(..)+                                        , Atom(..), Bond(..), GlobalID(..), LocalID(..)) +import           Control.Monad.IO.Class (MonadIO, liftIO)+import           Control.Exception      (evaluate)+import           Control.DeepSeq        (force, NFData) -oneModelSpecP :: Spec-oneModelSpecP = describe "One model." $-        it "correctly parses pdb with only one model without strings \"MODEL\" & \"ENDMDL\"" $ do-        let mt  = fromTextPDB . lenghtenLines $ T.pack ( "HEADER header\n" ++-                                         "TITLE STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME\n" ++-                                         "COMPND compnd\n" ++-                                         "SOURCE source\n" ++-                                         "KEYWDS keywds\n" ++-                                         "AUTHOR  Masha\n" ++-                                         "REVDAT revdat\n" ++-                                         "REMARK   1 REFERENCE 1\n" ++-                                         "SEQRES seqres\n" ++-                                         "CRYST1 cryst1\n" ++-                                         "ORIGX1 origx1 n=1\n" ++-                                         "SCALE2 sclaen n=2\n" ++-                                         "ATOM   2032  OXT CYS A 214      -4.546 -29.673  26.796  1.0 143.51           O  \n" ++-                                         "ATOM   2033  H   CYS A 214      -6.124 -27.225  26.558  1.00 15.00           H  \n" ++-                                         "TER    2034      CYS A 214                                                      \n" ++-                                         "ATOM   2035  N   GLU B   1      18.637 -61.583  66.852  1.0 118.48           N  \n" ++-                                         "TER   12534      ARG D 474                                                      \n" ++-                                         "HETATM12535  C1  NAG B 475       5.791 -20.194  -7.051  1.00 34.66           C  \n" ++-                                         "HETATM12538  C4  NAG B 475       6.943 -19.507  -9.597  1.00 25.87           C  \n" ++-                                         "ATOM   2036  CA  GLU B   1      19.722 -62.606  66.868  1.00 19.77           C  \n" ++-                                         "CONECT conect\n" ++-                                         "CONECT conect conect\n" ++-                                         "MASTER 1 2 3 4 5 6 7 8\n" ++-                                         "END"-                                       )-        mt `shouldBe` Right ([], Right oneModelPDB)+import           Data.Either            (fromRight)+import           Data.Vector            (Vector)+import qualified Data.Vector as V       (head, length, toList, concatMap)+import           Data.List              (find)+import           Data.Text              (Text)+import qualified Data.Text as T         (pack)+import           Data.Map.Strict        (Map, (!))+import qualified Data.Map.Strict as M   (fromList)+import           Data.Set               (Set)+import qualified Data.Set as S          (fromList, size, difference) -oneModelPDB :: PDB-oneModelPDB =  PDB { title = "STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME"-                   , remarks = Data.Map.Strict.singleton (Just 1) (V.singleton "REFERENCE 1")-                   , models = V.singleton $ V.fromList [V.fromList [Atom {atomSerial = 2032, atomName = " OXT", atomAltLoc = ' ', atomResName = "CYS", atomChainID = 'A', atomResSeq = 214, atomICode = ' ',-                                                                           atomX = -4.546, atomY = -29.673, atomZ = 26.796, atomOccupancy = 1.0, atomTempFactor = 143.51, atomElement = " O",-                                                                           atomCharge = "  "},-                                                                           Atom {atomSerial = 2033, atomName = " H  ", atomAltLoc = ' ', atomResName = "CYS", atomChainID = 'A',-                                                                           atomResSeq = 214, atomICode = ' ', atomX = -6.124, atomY = -27.225, atomZ = 26.558, atomOccupancy = 1.0, atomTempFactor = 15.0,-                                                                           atomElement = " H", atomCharge = "  "}],-                                                                           V.fromList [Atom {atomSerial = 2035, atomName = " N  ", atomAltLoc = ' ', atomResName = "GLU", atomChainID = 'B', atomResSeq = 1, atomICode = ' ',-                                                                           atomX = 18.637, atomY = -61.583, atomZ = 66.852, atomOccupancy = 1.0, atomTempFactor = 118.48, atomElement = " N", atomCharge = "  "},-                                                                           Atom {atomSerial = 2036, atomName = " CA ", atomAltLoc = ' ', atomResName = "GLU", atomChainID = 'B', atomResSeq = 1, atomICode = ' ',-                                                                           atomX = 19.722, atomY = -62.606, atomZ = 66.868, atomOccupancy = 1.0, atomTempFactor = 19.77, atomElement = " C", atomCharge = "  "}]-                                                                        ]-                  , otherFields = Data.Map.Strict.fromList [(HEADER, V.fromList [" header"]),(COMPND,V.fromList[" compnd"]),(SOURCE,V.fromList[" source"]),(KEYWDS,V.fromList[" keywds"]),(AUTHOR,V.fromList["  Masha"]),(REVDAT,V.fromList[" revdat"]),(SEQRES,V.fromList [" seqres"]),-                                                            (CRYST1,V.fromList [" cryst1"]),(ORIGX1,V.fromList[" origx1 n=1"]),(SCALE2,V.fromList [" sclaen n=2"]),(MASTER,V.fromList[" 1 2 3 4 5 6 7 8"])]-                   }+import           Test.Hspec -manyModelsSpecP :: Spec-manyModelsSpecP = describe "Some models." $-        it "correctly parses pdb with many models - they have strings \"MODEL\" & \"ENDMDL\" and text has disordered strings" $ do-        let mt = fromTextPDB . lenghtenLines $ T.pack ( "HEADER header\n" ++-                                        "TITLE STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME\n" ++-                                        "COMPND compnd\n" ++-                                        "SOURCE source\n" ++-                                        "REMARK 1   REFERENCE 1\n" ++-                                        "REMARK 2   Reference 2_1/1\n" ++-                                        "CRYST1 cryst1\n" ++-                                        "ORIGX1 origx1 n=1\n" ++-                                        "SCALE2 sclaen n=2\n" ++-                                        "MODEL  4  \n" ++-                                        "ATOM   2033  H   CYS A 214      -6.124 -27.225  26.558  1.00 15.00           H  \n" ++-                                        "TER    2034      CYS A 214                                                      \n" ++-                                        "ANISOU anisou\n" ++-                                        "ENDMDL \n" ++-                                        "MODEL  5 \n" ++-                                        "ATOM   2035  N   GLU B   1      18.637-691.583  66.852  1.0 118.48           N  \n" ++-                                        "ATOM  12531 HH12 ARG D 474      45.558 -39.551 -49.936  1.00 15.00           H  \n" ++-                                        "ENDMDL \n" ++-                                        "MODEL  6 \n" ++-                                        "ATOM  12532 HH21 ARG D 474      47.457 -38.007 -47.445  1.00 15.00           H  \n" ++-                                        "TER   12534      ARG D 474                                                      \n" ++-                                        "ATOM  12533 HH22 ARG D 474      47.405 -39.268 -48.629  1.00 15.00           H  \n" ++-                                        "HETATM12535  C1  NAG B 475       5.791 -20.194  -7.051  1.00 34.66           C  \n" ++-                                        "HETATM12538  C4  NAG B 475       6.943 -19.507  -9.597  1.00 25.87           C  \n" ++-                                        "CONECT conect conect\n" ++-                                        "CONECT conect\n" ++-                                        "ENDMDL\n" ++-                                        "SEQRES seqres\n" ++-                                        "KEYWDS keywds\n" ++-                                        "EXPDTA expdta\n" ++-                                        "AUTHOR  Masha\n" ++-                                        "REVDAT revdat\n" ++-                                        "MASTER 1 2 3 4 5 6 7 8\n"-                                      )-        mt `shouldBe` Right ([], Right manyModelsPDB)+rawPDBToModelConversionSingleChainSpec :: SpecWith ()+rawPDBToModelConversionSingleChainSpec = describe "Cobot Model from raw single chain PDB" $ do+  modelFromPDB <- runIO $ firstPDBModel "test/PDB/1PPE_I.pdb"+  modelFromMae <- runIO $ firstMaeModel "test/PDB/1PPE_I.mae"+  let (pdbBondCount, _, pdbChainCount, pdbAtomCount) = getStats modelFromPDB+  let (maeBondCount, _, _, _) = getStats modelFromMae -manyModelsPDB :: PDB-manyModelsPDB = PDB { title = "STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME"-                    , remarks = Data.Map.Strict.fromList [(Just 1, V.singleton "REFERENCE 1"), (Just 2, V.singleton "Reference 2_1/1")]-                    , models = V.fromList [V.fromList [V.fromList [Atom {atomSerial = 2033, atomName = " H  ", atomAltLoc = ' ', atomResName = "CYS",-                                                                                                  atomChainID = 'A', atomResSeq = 214, atomICode = ' ', atomX = -6.124, atomY = -27.225,-                                                                                                  atomZ = 26.558, atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}]-                                                                          ],-                                                    V.fromList [V.fromList [Atom {atomSerial = 2035, atomName = " N  ", atomAltLoc = ' ', atomResName = "GLU",-                                                                                                 atomChainID = 'B', atomResSeq = 1, atomICode = ' ', atomX = 18.637, atomY = -691.583,-                                                                                                 atomZ = 66.852, atomOccupancy = 1.0, atomTempFactor = 118.48, atomElement = " N", atomCharge = "  "}],-                                                                          V.fromList [Atom {atomSerial = 12531, atomName = "HH12", atomAltLoc = ' ', atomResName = "ARG",-                                                                                                 atomChainID = 'D', atomResSeq = 474, atomICode = ' ', atomX = 45.558, atomY = -39.551, atomZ = -49.936,-                                                                                                 atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}]-                                                                         ],-                                                    V.fromList [V.fromList [Atom {atomSerial = 12532, atomName = "HH21", atomAltLoc = ' ', atomResName = "ARG",-                                                                                                 atomChainID = 'D', atomResSeq = 474, atomICode = ' ', atomX = 47.457, atomY = -38.007, atomZ = -47.445,-                                                                                                 atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "},-                                                                                                Atom {atomSerial = 12533, atomName = "HH22", atomAltLoc = ' ', atomResName = "ARG",-                                                                                                 atomChainID = 'D', atomResSeq = 474, atomICode = ' ', atomX = 47.405, atomY = -39.268, atomZ = -48.629,-                                                                                                 atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}]-                                                                          ]-                                                      ]-                    , otherFields = Data.Map.Strict.fromList [(HEADER, V.fromList [" header"]),(COMPND,V.fromList[" compnd"]),(SOURCE,V.fromList[" source"]),(KEYWDS,V.fromList[" keywds"]),(EXPDTA,V.fromList[" expdta"]),(AUTHOR,V.fromList["  Masha"]),(REVDAT,V.fromList[" revdat"]),(SEQRES,V.fromList [" seqres"]),-                                                              (CRYST1,V.fromList [" cryst1"]),(ORIGX1,V.fromList[" origx1 n=1"]),(SCALE2,V.fromList [" sclaen n=2"]),(MASTER,V.fromList[" 1 2 3 4 5 6 7 8"])]-                    }+  it "Should have correct number of atoms" $ pdbAtomCount `shouldBe` 436+  it "Should have correct number of chains" $ pdbChainCount `shouldBe` 1+  it "Should restore bonds correctly" $ pdbBondCount `shouldBe` maeBondCount -noModelsSpecP :: Spec-noModelsSpecP = describe "No models." $-        it "correctly parses pdb without models (no ATOM, TER, HETATM strings)" $ do-        let mt = fromTextPDB . lenghtenLines $ T.pack ( "HEADER header\n" ++-                                        "TITLE STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME\n" ++-                                        "COMPND compnd\n" ++-                                        "SOURCE source\n" ++-                                        "KEYWDS keywds\n" ++-                                        "EXPDTA expdta\n" ++-                                        "AUTHOR  Masha\n" ++-                                        "REVDAT revdat\n" ++-                                        "REMARK 1   REFERENCE 1\n" ++-                                        "SEQRES seqres\n" ++-                                        "CRYST1 cryst1\n" ++-                                        "ORIGX1 origx1 n=1\n" ++-                                        "SCALE2 sclaen n=2\n" ++-                                        "MASTER 1 2 3 4 5 6 7 8\n"-                                      )-        mt `shouldBe` Right ([], Right noModelsPDB)+-- tripeptides are not checked as in checkBiggerMolecule +-- because there are inconsistencies in atom numbers between pdb and mae+bondsRestoringTripeptideSpec :: SpecWith ()+bondsRestoringTripeptideSpec = describe "Bonds should be restored correctly in tripeptides" $+  sequence_ $ checkTripeptide <$> tripeptides+    where+      tripeptides :: [String]+      tripeptides = ["ALA_3", "ARG_3", "ASN_3", "ASP_3", "CYS_3", "GLN_3", "GLU_3", "GLY_3", "HID_3", "HIE_3", "HIP_3", +                     "ILE_3", "LEU_3", "LYS_3", "MET_3", "PHE_3", "PRO_3", "SER_3", "THR_3", "TRP_3", "TYR_3", "VAL_3"]+      checkTripeptide :: String -> SpecWith (Arg Expectation)+      checkTripeptide tripeptideName = do+        modelFromMae <- runIO . firstMaeModel $ "test/PDB/BondsRestoring/" ++ tripeptideName ++ ".mae"+        modelFromPDB <- runIO . firstPDBModel $ "test/PDB/BondsRestoring/" ++ tripeptideName ++ ".pdb" -noModelsPDB :: PDB-noModelsPDB = PDB { title = "STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME"-                  , remarks = Data.Map.Strict.singleton (Just 1) (V.singleton "REFERENCE 1")-                  , models = V.empty-                  , otherFields = Data.Map.Strict.fromList [(HEADER, V.fromList [" header"]),(COMPND,V.fromList[" compnd"]),(SOURCE,V.fromList[" source"]),(KEYWDS,V.fromList[" keywds"]),(EXPDTA,V.fromList[" expdta"]),(AUTHOR,V.fromList["  Masha"]),(REVDAT,V.fromList[" revdat"]),(SEQRES,V.fromList [" seqres"]),-                                                            (CRYST1,V.fromList [" cryst1"]),(ORIGX1,V.fromList[" origx1 n=1"]),(SCALE2,V.fromList [" sclaen n=2"]),(MASTER,V.fromList[" 1 2 3 4 5 6 7 8"])]-                  }+        let pdbBondCount = V.length $ modelBonds modelFromPDB+        let maeBondCount = V.length $ modelBonds modelFromMae -allFieldsModelSpecP :: Spec-allFieldsModelSpecP = describe "PDB with all strings." $-        it "correctly parses pdb with all types of string" $ do-        let mt = fromTextPDB . lenghtenLines $ T.pack ( "HEADER header\n" ++-                                        "OBSLTE obslte\n" ++-                                        "TITLE STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME\n" ++-                                        "SPLIT split\n" ++-                                        "CAVEAT caveat\n" ++-                                        "COMPND compnd\n" ++-                                        "SOURCE source\n" ++-                                        "KEYWDS keywds\n" ++-                                        "EXPDTA expdta1\n" ++-                                        "EXPDTA expdta2\n" ++-                                        "NUMMDL nummdl\n" ++-                                        "MDLTYP mdltyp mdltyp\n" ++-                                        "AUTHOR  Masha\n" ++-                                        "REVDAT revdat\n" ++-                                        "SPRSDE sprsde\n" ++-                                        "JRNL   jrnl\n" ++-                                        "REMARK 1   REFERENCE 1\n" ++-                                        "DBREF dbref\n" ++-                                        "DBREF1 dbref1\n" ++-                                        "DBREF2 dbref2_1\n" ++-                                        "DBREF2 dbref2_2\n" ++-                                        "SEQADV seqadv\n" ++-                                        "SEQRES seqres\n" ++-                                        "MODRES modres\n" ++-                                        "HET    het\n" ++-                                        "HETNAM hetnam\n" ++-                                        "HETSYN hetsyn\n" ++-                                        "FORMUL of love\n" ++-                                        "HELIX helix\n" ++-                                        "SHEET sheet\n" ++-                                        "SSBOND ssbond\n" ++-                                        "LINK   link\n" ++-                                        "CISPEP cispep\n" ++-                                        "SITE   site\n" ++-                                        "CRYST1 cryst1\n" ++-                                        "ORIGX1 origx1 n=1\n" ++-                                        "SCALE2 sclaen n=2\n" ++-                                        "MTRIX3 matrixn n=3\n" ++-                                        "MODEL 1\n" ++-                                        "ATOM  12532 HH21 ARG D 474      47.457 -38.007 -47.445  1.00 15.00           H  \n" ++-                                        "TER   12534      ARG D 474                                                      \n" ++-                                        "HETATM12535  C1  NAG B 475       5.791 -20.194  -7.051  1.00 34.66           C  \n" ++-                                        "CONECT conect\n" ++-                                        "CONECT conect conect\n" ++-                                        "ENDMDL\n" ++-                                        "MASTER 1 2 3 4 5 6 7 8\n" ++-                                        "END"-                                      )-        mt `shouldBe` Right ([], Right pdbWithAllFields)+        it (tripeptideName ++ " equal bond count in Mae and PDB") $ pdbBondCount `shouldBe` maeBondCount -pdbWithAllFields :: PDB-pdbWithAllFields = PDB { title = "STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME"-                       , models = V.fromList [V.fromList [V.singleton Atom {atomSerial = 12532, atomName = "HH21", atomAltLoc = ' ', atomResName = "ARG",-                                     atomChainID = 'D', atomResSeq = 474, atomICode = ' ', atomX = 47.457, atomY = -38.007, atomZ = -47.445,-                                     atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}]]-                       , remarks = Data.Map.Strict.singleton (Just 1) (V.singleton "REFERENCE 1")-                       , otherFields = Data.Map.Strict.fromList [(HEADER, V.fromList [" header"]),(OBSLTE, V.fromList [" obslte"]), (SPLIT, V.fromList ["split"]), (CAVEAT, V.fromList [" caveat"]),-                                                                 (COMPND, V.fromList [" compnd"]),(SOURCE, V.fromList [" source"]),(KEYWDS, V.fromList [" keywds"]),-                                                                 (EXPDTA, V.fromList [" expdta1"," expdta2"]),(NUMMDL, V.fromList [" nummdl"]),(MDLTYP, V.fromList [" mdltyp mdltyp"]),-                                                                 (AUTHOR, V.fromList ["  Masha"]),(REVDAT, V.fromList [" revdat"]),(SPRSDE, V.fromList [" sprsde"]),-                                                                 (JRNL, V.fromList [" jrnl"]),(DBREF, V.fromList ["dbref"]), (DBREF1, V.fromList [" dbref1"]), (DBREF2, V.fromList [" dbref2_1"," dbref2_2"]),-                                                                 (SEQADV, V.fromList [" seqadv"]),(SEQRES, V.fromList [" seqres"]),(MODRES, V.fromList [" modres"]),-                                                                 (HET, V.fromList [" het"]), (HETNAM, V.fromList [" hetnam"]), (HETSYN, V.fromList [" hetsyn"]),(FORMUL, V.fromList [" of love"]),(HELIX, V.fromList ["helix"]),-                                                                 (SHEET, V.fromList ["sheet"]),(SSBOND, V.fromList [" ssbond"]),(LINK, V.fromList [" link"]),-                                                                 (CISPEP, V.fromList [" cispep"]),(SITE,V.fromList [" site"]),(CRYST1, V.fromList [" cryst1"]),(MTRIX3, V.fromList [" matrixn n=3"]),-                                                                 (ORIGX1, V.fromList [" origx1 n=1"]),(SCALE2, V.fromList [" sclaen n=2"]),(MASTER, V.fromList [" 1 2 3 4 5 6 7 8"])]-                       }+bondsRestoringBiggerMoleculesSpec :: SpecWith ()+bondsRestoringBiggerMoleculesSpec = describe "Bonds should be restored correctly in bigger molecules" $ do+  checkBiggerMolecule "3mxw_ab_b"+  checkBiggerMolecule "1vfb_ab_b"+  checkBiggerMolecule "4dn4_ag_b"+  where+    checkBiggerMolecule moleculeName = do+      modelFromPDB <- runIO . firstPDBModel $ "test/PDB/BondsRestoring/" ++ moleculeName ++ ".pdb"+      modelFromMae <- runIO . firstMaeModel $ "test/PDB/BondsRestoring/" ++ moleculeName ++ ".mae"+      let (pdbGlobalBondCount, pdbLocalBondCount, _, _) = getStats modelFromPDB+      let (maeGlobalBondCount, maeLocalBondCount, _, _) = getStats modelFromMae +      it (moleculeName ++ " equal global bond count in Mae and PDB") $ pdbGlobalBondCount `shouldBe` maeGlobalBondCount+      it (moleculeName ++ " equal local bond count in Mae and PDB") $ pdbLocalBondCount `shouldBe` maeLocalBondCount -emptySpecP :: Spec-emptySpecP = describe "empty PDB." $-        it "correctly parses empty pdb" $ do-        let mt = fromTextPDB ""-        mt `shouldBe` Right ([], Right emptyPdb)+      it (moleculeName ++ " no dublicate bonds") $ length (doubleBonds modelFromPDB) `shouldBe` 0 -emptyPdb :: PDB-emptyPdb = PDB { title = ""-               , models = V.empty-               , remarks = Data.Map.Strict.empty-               , otherFields = Data.Map.Strict.empty-               }+      let _globalBondSetPDB = globalBondSet modelFromPDB+      let _globalBondSetMae = globalBondSet modelFromMae+      let diffMaePDBGlobal = S.difference _globalBondSetMae _globalBondSetPDB+      let diffPDBMaeGlobal = S.difference _globalBondSetPDB _globalBondSetMae -trashBetweenModelsSpecP :: Spec-trashBetweenModelsSpecP = describe "PDB has trash." $-        it "correctly parses pdb with trash string between models and other field strings" $ do-        let mt = fromTextPDB  . lenghtenLines $ T.pack ( "trash strings 1\n" ++-                                        "HEADER header\n" ++-                                        "TITLE STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME\n" ++-                                        "REMARK 1   REFERENCE 1\n" ++-                                        "trash strings 2\n" ++-                                        "REMARK 2   Reference 2_1/1\n" ++-                                        "MODEL  4  \n" ++-                                        "ATOM   2033  H   CYS A 214      -6.124 -27.225  26.558  1.00 15.00           H  \n" ++-                                        "TER    2034      CYS A 214                                                      \n" ++-                                        "ANISOU anisou\n" ++-                                        "ENDMDL \n" ++-                                        "CRYST1 cryst1\n" ++-                                        "trash strings 3\n" ++-                                        "SEQRES seqres\n" ++-                                        "MODEL  5 \n" ++-                                        "ATOM   2035  N   GLU B   1      18.637-691.583  66.852  1.0 118.48           N  \n" ++-                                        "ATOM  12531 HH12 ARG D 474      45.558 -39.551 -49.936  1.00 15.00           H  \n" ++-                                        "ENDMDL \n" ++-                                        "trash strings 4\n" ++-                                        "MODEL  6 \n" ++-                                        "ATOM  12532 HH21 ARG D 474      47.457 -38.007 -47.445  1.00 15.00           H  \n" ++-                                        "TER   12534      ARG D 474                                                      \n" ++-                                        "ATOM  12533 HH22 ARG D 474      47.405 -39.268 -48.629  1.00 15.00           H  \n" ++-                                        "HETATM12535  C1  NAG B 475       5.791 -20.194  -7.051  1.00 34.66           C  \n" ++-                                        "HETATM12538  C4  NAG B 475       6.943 -19.507  -9.597  1.00 25.87           C  \n" ++-                                        "CONECT conect conect\n" ++-                                        "ENDMDL\n" ++-                                        "MASTER 1 2 3 4 5 6 7 8\n" ++-                                        "trash strings 5\n"-                                       )-        mt `shouldBe` Left "There are trash strings between model strings"+      it (moleculeName ++ " difference in Mae and PDB global bond sets") $ S.size diffMaePDBGlobal `shouldBe` 0+      it (moleculeName ++ " difference in PDB and Mae global bond sets") $ S.size diffPDBMaeGlobal `shouldBe` 0 -pdbWithoutTrash :: PDB-pdbWithoutTrash = PDB { title = "STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME"-                      , remarks = Data.Map.Strict.fromList [(Just 1, V.singleton "REFERENCE 1"), (Just 2, V.singleton "Reference 2_1/1")]-                      , models = V.fromList [V.fromList [V.fromList [Atom {atomSerial = 2033, atomName = " H  ", atomAltLoc = ' ', atomResName = "CYS",-                                                                                                    atomChainID = 'A', atomResSeq = 214, atomICode = ' ', atomX = -6.124, atomY = -27.225,-                                                                                                    atomZ = 26.558, atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}]-                                                                            ],-                                                      V.fromList [V.fromList [Atom {atomSerial = 2035, atomName = " N  ", atomAltLoc = ' ', atomResName = "GLU",-                                                                                                   atomChainID = 'B', atomResSeq = 1, atomICode = ' ', atomX = 18.637, atomY = -691.583,-                                                                                                   atomZ = 66.852, atomOccupancy = 1.0, atomTempFactor = 118.48, atomElement = " N", atomCharge = "  "}],-                                                                            V.fromList [Atom {atomSerial = 12531, atomName = "HH12", atomAltLoc = ' ', atomResName = "ARG",-                                                                                                   atomChainID = 'D', atomResSeq = 474, atomICode = ' ', atomX = 45.558, atomY = -39.551, atomZ = -49.936,-                                                                                                   atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}]-                                                                           ],-                                                      V.fromList [V.fromList [Atom {atomSerial = 12532, atomName = "HH21", atomAltLoc = ' ', atomResName = "ARG",-                                                                                                   atomChainID = 'D', atomResSeq = 474, atomICode = ' ', atomX = 47.457, atomY = -38.007, atomZ = -47.445,-                                                                                                   atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "},-                                                                                                  Atom {atomSerial = 12533, atomName = "HH22", atomAltLoc = ' ', atomResName = "ARG",-                                                                                                   atomChainID = 'D', atomResSeq = 474, atomICode = ' ', atomX = 47.405, atomY = -39.268, atomZ = -48.629,-                                                                                                   atomOccupancy = 1.0, atomTempFactor = 15.0, atomElement = " H", atomCharge = "  "}]-                                                                            ]-                                                        ]-                      , otherFields = Data.Map.Strict.fromList [(HEADER, V.fromList [" header"]),(SEQRES,V.fromList [" seqres"]),-                                                                (CRYST1,V.fromList [" cryst1"]),(MASTER,V.fromList[" 1 2 3 4 5 6 7 8"])]-                      }+      let _localBondSetPDB = localBondSet modelFromPDB+      let _localBondSetMae = localBondSet modelFromMae+      let diffMaePDBLocal = S.difference _localBondSetMae _localBondSetPDB+      let diffPDBMaeLocal = S.difference _localBondSetPDB _localBondSetMae -onlyOneModelSpecP :: Spec-onlyOneModelSpecP = describe "Only One model." $-        it "correctly parses pdb with only one model without other field/title/trash/remarks strings" $ do-        let mt = fromTextPDB . lenghtenLines $ T.pack ( "ATOM   2032  OXT CYS A 214      -4.546 -29.673  26.796  1.0 143.51           O  \n" ++-                                         "ATOM   2033  H   CYS A 214      -6.124 -27.225  26.558  1.00 15.00           H  \n" ++-                                         "TER    2034      CYS A 214                                                      \n" ++-                                         "ATOM   2035  N   GLU B   1      18.637 -61.583  66.852  1.0 118.48           N  \n" ++-                                         "TER   12534      ARG D 474                                                      \n" ++-                                         "HETATM12535  C1  NAG B 475       5.791 -20.194  -7.051  1.00 34.66           C  \n" ++-                                         "HETATM12538  C4  NAG B 475       6.943 -19.507  -9.597  1.00 25.87           C  \n" ++-                                         "ATOM   2036  CA  GLU B   1      19.722 -62.606  66.868  1.00 19.77           C  \n"-                                        )-        mt `shouldBe` Right ([], Right onlyOneModelPDB)+      it (moleculeName ++ " difference in Mae and PDB local bond sets") $ S.size diffMaePDBLocal `shouldBe` 0+      it (moleculeName ++ " difference in PDB and Mae local bond sets") $ S.size diffPDBMaeLocal `shouldBe` 0 -onlyOneModelPDB :: PDB-onlyOneModelPDB = PDB { title = ""-                      , remarks = Data.Map.Strict.empty-                      , models = V.singleton $ V.fromList [V.fromList [Atom {atomSerial = 2032, atomName = " OXT", atomAltLoc = ' ', atomResName = "CYS", atomChainID = 'A', atomResSeq = 214, atomICode = ' ',-                                                                              atomX = -4.546, atomY = -29.673, atomZ = 26.796, atomOccupancy = 1.0, atomTempFactor = 143.51, atomElement = " O",-                                                                              atomCharge = "  "},-                                                                              Atom {atomSerial = 2033, atomName = " H  ", atomAltLoc = ' ', atomResName = "CYS", atomChainID = 'A',-                                                                              atomResSeq = 214, atomICode = ' ', atomX = -6.124, atomY = -27.225, atomZ = 26.558, atomOccupancy = 1.0, atomTempFactor = 15.0,-                                                                              atomElement = " H", atomCharge = "  "}],-                                                                              V.fromList [Atom {atomSerial = 2035, atomName = " N  ", atomAltLoc = ' ', atomResName = "GLU", atomChainID = 'B', atomResSeq = 1, atomICode = ' ',-                                                                              atomX = 18.637, atomY = -61.583, atomZ = 66.852, atomOccupancy = 1.0, atomTempFactor = 118.48, atomElement = " N", atomCharge = "  "},-                                                                              Atom {atomSerial = 2036, atomName = " CA ", atomAltLoc = ' ', atomResName = "GLU", atomChainID = 'B', atomResSeq = 1, atomICode = ' ',-                                                                              atomX = 19.722, atomY = -62.606, atomZ = 66.868, atomOccupancy = 1.0, atomTempFactor = 19.77, atomElement = " C", atomCharge = "  "}]-                                                                              ]-                      , otherFields = Data.Map.Strict.empty-                      }+    localBondSet :: Model -> Set (Text, Int, Int, Int) -- (ChainID, ResidueNumber, LocalFrom, LocalTo)+    localBondSet Model{..} = S.fromList $ do+      Chain{..} <- V.toList modelChains+      Residue{..} <- V.toList chainResidues+      Bond (LocalID from) (LocalID to) _ <- V.toList resBonds -repeatedStringsSpecP :: Spec-repeatedStringsSpecP = describe "PDB with repeated other field strings." $-        it "correctly parses pdb with repeated other field strings (SOURCE)" $ do-        let mt = fromTextPDB . lenghtenLines $ T.pack ( "HEADER header\n" ++-                                        "TITLE STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME\n" ++-                                        "COMPND compnd\n" ++-                                        "EXPDTA expdta\n" ++-                                        "AUTHOR  Masha\n" ++-                                        "SOURCE source1\n" ++-                                        "KEYWDS keywds\n" ++-                                        "REVDAT revdat\n" ++-                                        "REMARK 1   REFERENCE 1\n" ++-                                        "SEQRES seqres\n" ++-                                        "CRYST1 cryst1\n" ++-                                        "ORIGX1 origx1 n=1\n" ++-                                        "SOURCE source2\n" ++-                                        "SCALE2 sclaen n=2\n" ++-                                        "MASTER 1 2 3 4 5 6 7 8\n" ++-                                        "SOURCE source3\n" ++-                                        "END"-                                      )-        mt `shouldBe` Right ([], Right repeatedStringsPDB)+      [(chainName, resNumber, from, to), (chainName, resNumber, to, from)] -repeatedStringsPDB :: PDB-repeatedStringsPDB = PDB { title = "STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME"-                         , remarks = Data.Map.Strict.singleton (Just 1) (V.singleton "REFERENCE 1")-                         , models = V.empty-                         , otherFields = Data.Map.Strict.fromList [(HEADER, V.fromList [" header"]),(COMPND,V.fromList[" compnd"]),(SOURCE,V.fromList[" source1", " source2", " source3"]),(KEYWDS,V.fromList[" keywds"]),(EXPDTA,V.fromList[" expdta"]),(AUTHOR,V.fromList["  Masha"]),(REVDAT,V.fromList[" revdat"]),(SEQRES,V.fromList [" seqres"]),-                                                                   (CRYST1,V.fromList [" cryst1"]),(ORIGX1,V.fromList[" origx1 n=1"]),(SCALE2,V.fromList [" sclaen n=2"]),(MASTER,V.fromList[" 1 2 3 4 5 6 7 8"])]-                         }+    doubleBonds :: Model -> [Bond GlobalID]+    doubleBonds Model{..} = doubleBonds' (V.toList modelBonds) []+      where+        doubleBonds' :: [Bond GlobalID] -> [Bond GlobalID] -> [Bond GlobalID]+        doubleBonds' [] acc = acc+        doubleBonds' (b:bs) acc = doubleBonds' bs . maybe acc (:acc) $ find (bondsEqual b) bs+        bondsEqual :: Bond GlobalID -> Bond GlobalID -> Bool+        bondsEqual b1 b2 = bondStart b1 == bondStart b2 && bondEnd b1 == bondEnd b2 -emptyRemarkSpecP :: Spec-emptyRemarkSpecP = describe "PDB with repeated remark strings without code." $-        it "correctly parses pdb with repeated remark strings without code" $ do-        let mt = fromTextPDB . lenghtenLines $ T.pack ( "REMARK 111 remark111_1/2\n" ++-                                        "COMPND compnd\n" ++-                                        "SOURCE source1\n" ++-                                        "KEYWDS keywds\n" ++-                                        "REMARK 111 remark111_2/2\n" ++-                                        "REVDAT revdat\n" ++-                                        "REMARK 2   remark2_1/1\n" ++-                                        "SOURCE source2\n" ++-                                        "REMARK     empty remark 1/2\n" ++-                                        "END    \n" ++-                                        "REMARK     empty remark 2/2\n" ++-                                        "SCALE2 sclaen n=2\n"-                                      )-        mt `shouldBe` Right ([], Right pdbWithEmptyRemarks)+    globalBondSet :: Model -> Set (Text,Text)+    globalBondSet Model{..} = bondSet getGlobalID (chainsAtomMap modelChains) modelBonds -pdbWithEmptyRemarks :: PDB-pdbWithEmptyRemarks = PDB { title = ""-                         , remarks = Data.Map.Strict.fromList [(Just 2, V.fromList ["remark2_1/1"]), (Just 111, V.fromList ["remark111_1/2", "remark111_2/2"]), (Nothing, V.fromList ["empty remark 1/2", "empty remark 2/2"])]-                         , models = V.empty-                         , otherFields = Data.Map.Strict.fromList [(COMPND,V.fromList[" compnd"]),(SOURCE,V.fromList[" source1", " source2"]),-                                                                    (KEYWDS,V.fromList[" keywds"]), (REVDAT,V.fromList[" revdat"]), (SCALE2,V.fromList [" sclaen n=2"])]-                         }+    chainsAtomMap :: Vector Chain -> Map Int (Text, Atom)+    chainsAtomMap chains = M.fromList $ concatMap chainAtomPreMap chains -emptyModelSpecP :: Spec-emptyModelSpecP = describe "PDB with one empty model." $-        it "correctly parses pdb with one model without strings inside" $ do-        let mt = fromTextPDB . lenghtenLines $ T.pack ( "trash strings 1\n" ++-                                        "HEADER header\n" ++-                                        "TITLE STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME\n" ++-                                        "REMARK 1   REFERENCE 1\n" ++-                                        "trash strings 2\n" ++-                                        "REMARK 2   Reference 2_1/1\n" ++-                                        "MODEL  4  \n" ++-                                        "ENDMDL\n" ++-                                        "MASTER 1 2 3 4 5 6 7 8\n" ++-                                        "trash strings 5\n"-                                      )-        mt `shouldBe` Right ([], Right pdbEmptyModel)+    chainAtomPreMap :: Chain -> [(Int, (Text, Atom))]+    chainAtomPreMap Chain{..} = V.toList . fmap (\a -> (getGlobalID $ atomId a, (chainName, a))) $ V.concatMap resAtoms chainResidues +    bondSet :: (a -> Int) -> Map Int (Text, Atom) -> Vector (Bond a) -> Set (Text,Text)+    bondSet getID atomMap bonds = S.fromList $ do+      Bond{..} <- V.toList bonds+      let atomFromId = formAtomId $ atomMap ! getID bondStart+      let atomToId = formAtomId $ atomMap ! getID bondEnd+      [(atomFromId, atomToId), (atomToId, atomFromId)]+      where+        formAtomId :: (Text, Atom) -> Text+        formAtomId (chainId, Atom{..}) = chainId <> "_" <> atomName <> "_" <> T.pack (show $ getGlobalID atomId) -pdbEmptyModel :: PDB-pdbEmptyModel = PDB { title = "STRUCTURE OF THE TRANSFORMED MONOCLINIC  LYSOZYME"-                    , remarks = Data.Map.Strict.fromList [(Just 1, V.singleton "REFERENCE 1"), (Just 2, V.singleton "Reference 2_1/1")]-                    , models = V.empty-                    , otherFields = Data.Map.Strict.fromList [(HEADER, V.fromList [" header"]),-                                                              (MASTER,V.fromList[" 1 2 3 4 5 6 7 8"])]-                    } -lenghtenLines :: Text -> Text-lenghtenLines text = longLinedText-   where-       textLines = T.lines text-       longTextLines = changeLine <$> textLines-       desiredLength = 80  -- cause it is max length in standart pdb file-       longLinedText = T.intercalate "\n" longTextLines+getStats :: Model -> (Int, Int, Int, Int)+getStats model = (globalBondCount, localBondCount, chainCount, atomCount)+  where +    globalBondCount = V.length $ modelBonds model+    localBondCount = V.length . V.concatMap (V.concatMap resBonds . chainResidues) $ modelChains model+    chainCount = V.length $ modelChains model+    atomCount = V.length . V.concatMap (V.concatMap resAtoms . chainResidues) $ modelChains model -       changeLine :: Text -> Text-       changeLine line | T.length line > desiredLength = T.take desiredLength line-                       | T.length line < desiredLength = line <> T.replicate spacesCount " "-                       | otherwise = line-            where-                    spacesCount = desiredLength - T.length line+ef :: NFData a => a -> IO a+ef = evaluate . force +firstPDBModel :: (MonadIO m) => FilePath -> m Model+firstPDBModel filepath = do+  eitherPDB <- modelsFromPDBFile filepath+  let (_, models) = fromRight undefined eitherPDB+  liftIO . ef $ V.head models++firstMaeModel :: (MonadIO m) => FilePath -> m Model+firstMaeModel filepath = do+  eitherMae <- modelsFromMaeFile filepath+  -- `evaluate . force` fails for some reason+  pure . V.head $ fromRight undefined eitherMae

@@ -13,6 +13,7 @@ import           System.IO import           Test.Hspec import           UniprotSpec+import           PDBParserSpec  main :: IO () main = do@@ -53,5 +54,8 @@          repeatedStringsSpecP          emptyRemarkSpecP          emptyModelSpecP+         rawPDBToModelConversionSingleChainSpec+         bondsRestoringTripeptideSpec+         bondsRestoringBiggerMoleculesSpec          -- Structure          structureSpec