Mercurial > repos > yufei-luo > s_mart
diff commons/core/parsing/CoordsParser.py @ 6:769e306b7933
Change the repository level.
author | yufei-luo |
---|---|
date | Fri, 18 Jan 2013 04:54:14 -0500 |
parents | |
children |
line wrap: on
line diff
--- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/commons/core/parsing/CoordsParser.py Fri Jan 18 04:54:14 2013 -0500 @@ -0,0 +1,137 @@ +# +# Copyright INRA-URGI 2009-2010 +# +# This software is governed by the CeCILL license under French law and +# abiding by the rules of distribution of free software. You can use, +# modify and/ or redistribute the software under the terms of the CeCILL +# license as circulated by CEA, CNRS and INRIA at the following URL +# "http://www.cecill.info". +# +# As a counterpart to the access to the source code and rights to copy, +# modify and redistribute granted by the license, users are provided only +# with a limited warranty and the software's author, the holder of the +# economic rights, and the successive licensors have only limited +# liability. +# +# In this respect, the user's attention is drawn to the risks associated +# with loading, using, modifying and/or developing or reproducing the +# software by the user in light of its specific status of free software, +# that may mean that it is complicated to manipulate, and that also +# therefore means that it is reserved for developers and experienced +# professionals having in-depth computer knowledge. Users are therefore +# encouraged to load and test the software's suitability as regards their +# requirements in conditions enabling the security of their systems and/or +# data to be ensured and, more generally, to use and operate it in the +# same conditions as regards security. +# +# The fact that you are presently reading this means that you have had +# knowledge of the CeCILL license and that you accept its terms. +# +import re +import sys +from SMART.Java.Python.structure.Mapping import Mapping +from commons.core.parsing.MapperParser import MapperParser +from SMART.Java.Python.structure.SubMapping import SubMapping +from SMART.Java.Python.misc import Utils + +class CoordsParser(MapperParser): + """A class that parses the .coords output of Nucmer""" + + def __init__(self, fileName, verbosity = 0): + self._lineParseRe = re.compile(r"^\s*(?P<tStart>\d+)\s+(?P<tEnd>\d+)\s+\|\s+(?P<qStart>\d+)\s+(?P<qEnd>\d+)\s+\|\s+(?P<tLength>\d+)\s+(?P<qLength>\d+)\s+\|\s+(?P<identity>\d+\.?\d*)\s+\|\s+(?P<tName>[\w\|\:\-]+)\s+(?P<qName>.*)\s*$") + self._lineParseRe2 = re.compile(r"^\s*(?P<tStart>\d+)\s+(?P<tEnd>\d+)\s+(?P<qStart>\d+)\s+(?P<qEnd>\d+)\s+(?P<tLength>\d+)\s+(?P<qLength>\d+)\s+(?P<identity>\d+\.?\d*)\s+(?P<rlen>\d+\.?\d*)\s+(?P<qlen>\d+\.?\d*)\s+(?P<rcov>\d+\.?\d*)\s+(?P<qcov>\d+\.?\d*)\s+(?P<rframe>[-]?\d+\.?\d*)\s+(?P<qframe>[-]?\d+\.?\d*)\s+(?P<tName>[\w\|\:\-]+)\s+(?P<qName>.*)\s*$") + self._lineParseRe3 = re.compile(r"^\s*(?P<tStart>\d+)\s+(?P<tEnd>\d+)\s+\|\s+(?P<qStart>\d+)\s+(?P<qEnd>\d+)\s+\|\s+(?P<tLength>\d+)\s+(?P<qLength>\d+)\s+\|\s+(?P<identity>\d+\.?\d*)\s+(?P<sim>\d+\.?\d*)\s+(?P<stp>\d+\.?\d*)\s+\|\s+(?P<rframe>[-]?\d+\.?\d*)\s+(?P<qframe>[-]?\d+\.?\d*)\s+(?P<tName>[\w\|\:\-]+)\s+(?P<qName>.*)\s*$") + self._lineParseRe4 = re.compile(r"^\s*(?P<tStart>\d+)\s+(?P<tEnd>\d+)\s+(?P<qStart>\d+)\s+(?P<qEnd>\d+)\s+(?P<tLength>\d+)\s+(?P<qLength>\d+)\s+(?P<identity>\d+\.?\d*)\s+(?P<sim>\d+\.?\d*)\s+(?P<stp>\d+\.?\d*)\s+(?P<rlen>\d+\.?\d*)\s+(?P<qlen>\d+\.?\d*)\s+(?P<rcov>\d+\.?\d*)\s+(?P<qcov>\d+\.?\d*)\s+(?P<rframe>[-]?\d+\.?\d*)\s+(?P<qframe>[-]?\d+\.?\d*)\s+(?P<tName>[\w\|\:\-]+)\s+(?P<qName>.*)\s*$") + self.lineType = 1 + MapperParser.__init__(self, fileName, verbosity) + + def getFileFormats(): + return ["coords"] + getFileFormats = staticmethod(getFileFormats) + + def skipFirstLines(self): + while True: + line = self.handle.readline() + self.currentLineNb += 1 + if line == "": + break + if "=====" in line: + break + if "[S1]\t[E1]\t[S2]\t[E2]\t[LEN 1]\t[LEN 2]\t[% IDY]\t[LEN R]\t[LEN Q]\t[COV R]\t[COV Q]\t[FRM]\t[TAGS]" in line: + self.lineType = 2 + break + if "[S1] [E1] | [S2] [E2] | [LEN 1] [LEN 2] | [% IDY] [% SIM] [% STP] | [FRM] [TAGS]" in line: + self.lineType = 3 + + if "[% IDY]\t[% SIM]\t[% STP]" in line and "[LEN Q]"in line: + self.lineType = 4 + break + + def parseLine(self, line): + + if self.lineType == 1 : + m = self._lineParseRe.search(line) + elif self.lineType == 2: + m = self._lineParseRe2.search(line) + elif self.lineType == 3: + m = self._lineParseRe3.search(line) + elif self.lineType == 4: + m = self._lineParseRe4.search(line) + if m == None: + sys.exit("\nLine %d '%s' does not have a NucMer format" % (self.currentLineNb, line)) + + mapping = Mapping() + + subMapping = SubMapping() + subMapping.queryInterval.setName(m.group("qName")) + subMapping.queryInterval.setStart(min(int(m.group("qStart")), int(m.group("qEnd")))) + subMapping.queryInterval.setEnd(max(int(m.group("qStart")), int(m.group("qEnd")))) + subMapping.queryInterval.setSize(int(m.group("qLength"))) + subMapping.queryInterval.setDirection(int(m.group("qEnd")) - int(m.group("qStart"))) + + subMapping.targetInterval.setChromosome(m.group("tName")) + subMapping.targetInterval.setStart(min(int(m.group("tStart")), int(m.group("tEnd")))) + subMapping.targetInterval.setEnd(max(int(m.group("tStart")), int(m.group("tEnd")))) + subMapping.targetInterval.setSize(int(m.group("tLength"))) + subMapping.targetInterval.setDirection(int(m.group("tEnd")) - int(m.group("tStart"))) + + subMapping.setDirection(int(m.group("qEnd")) - int(m.group("qStart"))) + subMapping.setSize(min(int(m.group("qLength")), int(m.group("tLength")))) + subMapping.setIdentity(float(m.group("identity"))) + + mapping.addSubMapping(subMapping) + mapping.targetInterval.setStart(min(int(m.group("tStart")), int(m.group("tEnd")))) + mapping.targetInterval.setEnd(max(int(m.group("tStart")), int(m.group("tEnd")))) + mapping.targetInterval.setSize(int(m.group("tLength"))) + mapping.targetInterval.setChromosome(m.group("tName")) + + mapping.queryInterval.setStart(min(int(m.group("qStart")), int(m.group("qEnd")))) + mapping.queryInterval.setEnd(max(int(m.group("qStart")), int(m.group("qEnd")))) + mapping.queryInterval.setSize(int(m.group("qLength"))) + mapping.queryInterval.setName(m.group("qName")) + mapping.setDirection(int(m.group("qEnd")) - int(m.group("qStart"))) + mapping.setSize(min(int(m.group("qLength")), int(m.group("tLength")))) + mapping.setIdentity(float(m.group("identity"))) + mapping.setTagValue("feature", "match") + mapping.setTagValue("Target", "%s %d %d" % (m.group("qName"), int(m.group("qStart")), int(m.group("qEnd")))) + + if self.lineType ==2 or self.lineType ==4: + mapping.setTagValue("target_pident", float(m.group("identity"))) + mapping.setTagValue("target_pcover", float(m.group("qcov"))) + mapping.setTagValue("target_length", int(m.group("qlen"))) + + +# Specific to Mark Work. Commented lines because of possible slowdown. +# for line in self.handle: +# string1 = line.strip() +# self.currentLineNb += 1 +# break +# for line in self.handle: +# string2 = line.strip() +# self.currentLineNb += 1 +# break +# print(len(string1),len(string2)) +# mapping.setNbMismatches(Utils.getHammingDistance(string1, string2)) + mapping.setNbGaps(0) + + return mapping