Back to index

python-biopython  1.60
Public Member Functions | Public Attributes | Private Attributes
Bio.SeqIO.QualityIO.QualPhredWriter Class Reference
Inheritance diagram for Bio.SeqIO.QualityIO.QualPhredWriter:
Inheritance graph
[legend]
Collaboration diagram for Bio.SeqIO.QualityIO.QualPhredWriter:
Collaboration graph
[legend]

List of all members.

Public Member Functions

def __init__
def write_record
def write_header
def write_footer
def write_records
def write_file
def clean

Public Attributes

 wrap
 record2title
 handle

Private Attributes

 _record_written

Detailed Description

Class to write QUAL format files (using PHRED quality scores).

Although you can use this class directly, you are strongly encouraged
to use the Bio.SeqIO.write() function instead.  For example, this code
reads in a FASTQ file and saves the quality scores into a QUAL file:

>>> from Bio import SeqIO
>>> record_iterator = SeqIO.parse(open("Quality/example.fastq"), "fastq")
>>> out_handle = open("Quality/temp.qual", "w")
>>> SeqIO.write(record_iterator, out_handle, "qual")
3
>>> out_handle.close()

This code is also called if you use the .format("qual") method of a
SeqRecord.

P.S. Don't forget to clean up the temp file if you don't need it anymore:

>>> import os
>>> os.remove("Quality/temp.qual")

Definition at line 1443 of file QualityIO.py.


Constructor & Destructor Documentation

def Bio.SeqIO.QualityIO.QualPhredWriter.__init__ (   self,
  handle,
  wrap = 60,
  record2title = None 
)
Create a QUAL writer.

Arguments:
 - handle - Handle to an output file, e.g. as returned
    by open(filename, "w")
 - wrap   - Optional line length used to wrap sequence lines.
    Defaults to wrapping the sequence at 60 characters
    Use zero (or None) for no wrapping, giving a single
    long line for the sequence.
 - record2title - Optional function to return the text to be
    used for the title line of each record.  By default
    a combination of the record.id and record.description
    is used.  If the record.description starts with the
    record.id, then just the record.description is used.

The record2title argument is present for consistency with the
Bio.SeqIO.FastaIO writer class.

Definition at line 1465 of file QualityIO.py.

01465 
01466     def __init__(self, handle, wrap=60, record2title=None):
01467         """Create a QUAL writer.
01468 
01469         Arguments:
01470          - handle - Handle to an output file, e.g. as returned
01471                     by open(filename, "w")
01472          - wrap   - Optional line length used to wrap sequence lines.
01473                     Defaults to wrapping the sequence at 60 characters
01474                     Use zero (or None) for no wrapping, giving a single
01475                     long line for the sequence.
01476          - record2title - Optional function to return the text to be
01477                     used for the title line of each record.  By default
01478                     a combination of the record.id and record.description
01479                     is used.  If the record.description starts with the
01480                     record.id, then just the record.description is used.
01481 
01482         The record2title argument is present for consistency with the
01483         Bio.SeqIO.FastaIO writer class.
01484         """
01485         SequentialSequenceWriter.__init__(self, handle)
01486         #self.handle = handle
01487         self.wrap = None
01488         if wrap:
01489             if wrap < 1:
01490                 raise ValueError
01491         self.wrap = wrap
01492         self.record2title = record2title


Member Function Documentation

def Bio.SeqIO.Interfaces.SequenceWriter.clean (   self,
  text 
) [inherited]
Use this to avoid getting newlines in the output.

Definition at line 166 of file Interfaces.py.

00166 
00167     def clean(self, text):
00168         """Use this to avoid getting newlines in the output."""
00169         return text.replace("\n", " ").replace("\r", " ").replace("  ", " ")
    

Here is the caller graph for this function:

def Bio.SeqIO.Interfaces.SequentialSequenceWriter.write_file (   self,
  records 
) [inherited]
Use this to write an entire file containing the given records.

records - A list or iterator returning SeqRecord objects

This method can only be called once.  Returns the number of records
written.

Reimplemented from Bio.SeqIO.Interfaces.SequenceWriter.

Definition at line 262 of file Interfaces.py.

00262 
00263     def write_file(self, records):
00264         """Use this to write an entire file containing the given records.
00265 
00266         records - A list or iterator returning SeqRecord objects
00267 
00268         This method can only be called once.  Returns the number of records
00269         written.
00270         """
00271         self.write_header()
00272         count = self.write_records(records)
00273         self.write_footer()
00274         return count

Here is the call graph for this function:

Definition at line 218 of file Interfaces.py.

00218 
00219     def write_footer(self):
00220         assert self._header_written, "You must call write_header() first"
00221         assert self._record_written, "You have not called write_record() or write_records() yet"
00222         assert not self._footer_written, "You have aleady called write_footer()"
00223         self._footer_written = True

Here is the caller graph for this function:

Definition at line 212 of file Interfaces.py.

00212 
00213     def write_header(self):
00214         assert not self._header_written, "You have aleady called write_header()"
00215         assert not self._record_written, "You have aleady called write_record() or write_records()"
00216         assert not self._footer_written, "You have aleady called write_footer()"
00217         self._header_written = True
        

Here is the caller graph for this function:

Write a single QUAL record to the file.

Reimplemented from Bio.SeqIO.Interfaces.SequentialSequenceWriter.

Definition at line 1493 of file QualityIO.py.

01493 
01494     def write_record(self, record):
01495         """Write a single QUAL record to the file."""
01496         assert self._header_written
01497         assert not self._footer_written
01498         self._record_written = True
01499 
01500         handle = self.handle
01501         wrap = self.wrap
01502 
01503         if self.record2title:
01504             title = self.clean(self.record2title(record))
01505         else:
01506             id = self.clean(record.id)
01507             description = self.clean(record.description)
01508             if description and description.split(None, 1)[0]==id:
01509                 #The description includes the id at the start
01510                 title = description
01511             elif description:
01512                 title = "%s %s" % (id, description)
01513             else:
01514                 title = id
01515         handle.write(">%s\n" % title)
01516 
01517         qualities = _get_phred_quality(record)
01518         try:
01519             #This rounds to the nearest integer.
01520             #TODO - can we record a float in a qual file?
01521             qualities_strs = [("%i" % round(q, 0)) for q in qualities]
01522         except TypeError, e:
01523             if None in qualities:
01524                 raise TypeError("A quality value of None was found")
01525             else:
01526                 raise e
01527 
01528         if wrap > 5:
01529             #Fast wrapping
01530             data = " ".join(qualities_strs)
01531             while True:
01532                 if len(data) <= wrap:
01533                     self.handle.write(data + "\n")
01534                     break
01535                 else:
01536                     #By construction there must be spaces in the first X chars
01537                     #(unless we have X digit or higher quality scores!)
01538                     i = data.rfind(" ", 0, wrap)
01539                     handle.write(data[:i] + "\n")
01540                     data = data[i+1:]
01541         elif wrap:
01542             #Safe wrapping
01543             while qualities_strs:
01544                 line = qualities_strs.pop(0)
01545                 while qualities_strs \
01546                 and len(line) + 1 + len(qualities_strs[0]) < wrap:
01547                     line += " " + qualities_strs.pop(0)
01548                 handle.write(line + "\n")
01549         else:
01550             #No wrapping
01551             data = " ".join(qualities_strs)
01552             handle.write(data + "\n")

def Bio.SeqIO.Interfaces.SequentialSequenceWriter.write_records (   self,
  records 
) [inherited]
Write multiple record to the output file.

records - A list or iterator returning SeqRecord objects

Once you have called write_header() you can call write_record()
and/or write_records() as many times as needed.  Then call
write_footer() and close().

Returns the number of records written.

Definition at line 240 of file Interfaces.py.

00240 
00241     def write_records(self, records):
00242         """Write multiple record to the output file.
00243 
00244         records - A list or iterator returning SeqRecord objects
00245 
00246         Once you have called write_header() you can call write_record()
00247         and/or write_records() as many times as needed.  Then call
00248         write_footer() and close().
00249 
00250         Returns the number of records written.
00251         """
00252         #Default implementation:
00253         assert self._header_written, "You must call write_header() first"
00254         assert not self._footer_written, "You have already called write_footer()"
00255         count = 0
00256         for record in records:
00257             self.write_record(record)
00258             count += 1
00259         #Mark as true, even if there where no records
00260         self._record_written = True
00261         return count

Here is the call graph for this function:

Here is the caller graph for this function:


Member Data Documentation

Reimplemented from Bio.SeqIO.Interfaces.SequentialSequenceWriter.

Definition at line 1497 of file QualityIO.py.

Reimplemented from Bio.SeqIO.Interfaces.SequenceWriter.

Definition at line 207 of file Interfaces.py.

Definition at line 1491 of file QualityIO.py.

Definition at line 1486 of file QualityIO.py.


The documentation for this class was generated from the following file: