Source code for plot

#!/usr/bin/env python3.2
# -*- coding: utf-8 -*-
# encoding: utf-8
"""
Plotting functions for pymzML
"""

# pymzml
#
# Copyright (C) 2010-2011 T. Bald, J. Barth, A. Niehues, C. Fufezan
#
#    This program is free software: you can redistribute it and/or modify
#    it under the terms of the GNU General Public License as published by
#    the Free Software Foundation, either version 3 of the License, or
#    (at your option) any later version.
#
#    This program is distributed in the hope that it will be useful,
#    but WITHOUT ANY WARRANTY; without even the implied warranty of
#    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
#    GNU General Public License for more details.
#
#    You should have received a copy of the GNU General Public License
#    along with this program.  If not, see <http://www.gnu.org/licenses/>.

from __future__ import print_function
import math
import sys
import collections


[docs]class Factory(object): """ Class to plot pymzml.spec.Spectrum as svg/xhtml. :param filename: Name for the output file. Default = "spectra.xhtml" :type filename: string Example: >>> import pymzml, get_example_file >>> mzMLFile = 'profile-mass-spectrum.mzml' >>> example_file = get_example_file.open_example(mzMLFile) >>> run = pymzml.run.Run("../mzML_example_files/"+mzMLFile, precisionMSn = 250e-6) >>> p = pymzml.plot.Factory() >>> for spec in run: >>> p.newPlot() >>> p.add(spec.peaks, color=(200,00,00), style='circles') >>> p.add(spec.centroidedPeaks, color=(00,00,00), style='sticks') >>> p.add(spec.reprofiledPeaks, color=(00,255,00), style='circles') >>> p.save( filename="output/plotAspect.xhtml" , mzRange = [745.2,745.6] ) """ def __init__(self,filename = None): self.filename = filename if filename != None else "spectra.xhtml" self.plots = [ ] # list of plots, where each plot hold a list of mz-i lists that should be plotted in the same plot self.styles = [ ] self.colors = [ ] self.header = [ ] self.labelPos = [ ] self.mzRanges = [ ] self.normalizations = [] pass
[docs] def newPlot(self, header = None , mzRange = [None,None] , normalize = False): """ Add new plot to the plotFactory. :param header: an optional title that will be printed above the plot :type header: string :param mzRange: Boundaries of the new plot :type mzRange: tuple of minMZ,maxMZ :param normalize: whether or not the individal data sets are normalized in the plot :type boolean: """ self.plots.append( [] ) self.styles.append( [] ) self.colors.append( [] ) self.labelPos.append( set() ) self.mzRanges.append( mzRange ) self.normalizations.append( normalize ) if header == None: self.header.append( '' ) else: self.header.append( header ) return
[docs] def add(self,data, color=(00,00,00) , style = 'sticks', mzRange = [None,None] ): """ Add data to the graph. :param data: The data added to the graph :type data: list of tuples (mz,i) :param color: color encoded in RGB. Default = (0,0,0) :type color: tuple (R,G,B) :param style: plotting style. Default = "circles". :type style: string :param mzRange: Boundaries that should be added to the current plot :type mzRange: tuple of minMZ,maxMZ Currently supported styles are: * 'emptycircles' * 'circles' * 'sticks' * 'squares' * 'area' * 'triangle' * 'label' NOTE: The data format for label style is [( mz1, 'label1' ), ( mz2, 'label2' ), ( mz3, 'label3' ) ]. """ if len(self.plots) == 0: self.newPlot() self.styles[-1].append(style) self.colors[-1].append(color) selective = True if mzRange == [None,None]: if self.mzRanges[-1] == [None,None]: selective = False self.plots[-1].append(data) else: mzRange = self.mzRanges[-1] else: self.mzRanges[-1] = mzRange # NOTE this allows us to define the mzRange in the call of the add function, otherwise it has to be in the call of the new function # NOTE Problem can be: mzRanges[-1] can be overwritten ... don't know what is best here if selective: self.plots[-1].append([(mz,i) for mz,i in data if mzRange[0] <= mz <=mzRange[1]]) #print("Added {0} data points".format(len(self.plots[-1][-1]))) return
[docs] def info(self): """ Returns summary about the plotting factory, i.e.how many plots and how many datasets per plot. """ print(""" Factory holds {0} unique plots""".format(len(self.plots))) for i,plot in enumerate(self.plots): print(" Plot {0} holds {1} unique datasets".format(i,len(plot))) for j,dataset in enumerate(plot): print(" Dataset {0} holds {1} datapoints".format(j,len(dataset))) print() return
[docs] def save(self, filename = None, mzRange = [None,None]): """ Saves all plots and their data points that have been added to the plotFactory. :param filename: Name for the output file. Default = "spectra.xhtml" :type filename: string :param mzRange: m/z range which should be considered [start, end]. Default = [ ``None`` , ``None`` ] :type mzRange: list """ if filename == None: filename = self.filename io = open("{0}".format(filename), 'w' ) # NOTE: we need to check if file exists and/or if file can be written ... print("""<html xmlns="http://www.w3.org/1999/xhtml"> <head> <title>pymzml spectrum visualisation</title> <style> </style> </head> <body style="position:relative; z-index:0; width:100%; height:100%;"> """, file = io) for plotNumber,plot in enumerate(self.plots): w = 1000 h = 700 heada = self.header[plotNumber] resolved_mzRange = [ None, None ] if mzRange != [ None, None ]: resolved_mzRange = mzRange else: if self.mzRanges[plotNumber] == [None,None]: # we have to determine ranges ... for dataset in plot: print("Determing mzRanges ...") if len(dataset) == 0: continue if resolved_mzRange[0] == None: resolved_mzRange[0] = min(dataset)[0] else: if resolved_mzRange[0] > min(dataset)[0]: resolved_mzRange[0] = min(dataset)[0] if resolved_mzRange[1] == None: resolved_mzRange[1] = max(dataset)[0] else: if resolved_mzRange[1] < max(dataset)[0]: resolved_mzRange[1] = max(dataset)[0] else: # Ranges have been defined by newPlot() # i.e we want to plot no matter what ... resolved_mzRange = self.mzRanges[plotNumber] if resolved_mzRange == [None,None]: print("Skipping empty plot, plotNumber {0} , header :{1} ...".format(plotNumber,heada)) continue i_max = 0 i_min = 0 for dataset in plot: for mz,i in dataset: if resolved_mzRange[0] <= mz <= resolved_mzRange[1]: try: i = float(i) except: continue if i >= i_max: i_max = i if i_max == 0: i_max = 0.1 if self.normalizations[plotNumber] == True: resolved_i_max = 1.05 else: resolved_i_max = i_max * 1.05 print(""" <svg xmlns="http://www.w3.org/2000/svg" version="1.1" width="{w}" height="{h}" style="position:relative; top:0; left:0; z-index:-1;"> """.format(w = w, h = h), file = io) PADDING = (300,50,100,50) # css style convention baseline = h-PADDING[2] # Drawing Frame ... for x1,y1,x2,y2 in [ (PADDING[3],PADDING[0],w-PADDING[1],PADDING[0]), (w-PADDING[1],PADDING[0],w-PADDING[1],h-PADDING[2]), (PADDING[3],h-PADDING[2],w-PADDING[1],h-PADDING[2]), (PADDING[3],PADDING[0],PADDING[3],h-PADDING[2]), ]: print("""<line x1="{0}" y1="{1}" x2="{2}" y2="{3}" style="stroke:rgb(00,00,00);stroke-width:1"/>""".format(x1,y1,x2,y2),file=io) if heada != '': heada_y = 30 for subheada in heada.split("\n"): print("""<text x="{x}" y="{y}" font-family="Courier" text-anchor="start" font-size="{fontsize}" fill="black" >{0}</text> """.format(subheada, x = 0, y = heada_y, fontsize = 24), file = io) heada_y += 30 pixelpermz = float( w - PADDING[1] - PADDING[3] ) / float(resolved_mzRange[1]-resolved_mzRange[0]) # NOTE this returns ZeroDivisionError if the range to plot contains only one peak pixelper_i = float( h - PADDING[0] - PADDING[2] ) / float(resolved_i_max-i_min) delta = resolved_mzRange[1]-resolved_mzRange[0] sep = round(math.log(delta,10)) factor = 10**(-1*(sep-1)) # grid # Drawing ticks and x labels for _ in range(int(round(resolved_mzRange[0]*factor)),int(round(resolved_mzRange[1]*factor))+1 ): x = int(round(PADDING[3] + (float(_)/float(factor) - resolved_mzRange[0] ) * pixelpermz )) y = baseline color = 77 #if _ % 10 == 0 else 125 # Dashed lines # print("""<line x1="{0}" y1="{1}" x2="{0}" y2="{2}" # style="stroke-dasharray: 9, 5;stroke:rgb({3},{3},{3});stroke-width:0.5"/>""".format(x,PADDING[0],baseline,color),file=io) print("""<line x1="{0}" y1="{1}" x2="{0}" y2="{2}" style="stroke:rgb(00,00,00);stroke-width:0.5;"/>""".format(x,baseline,baseline+10),file=io) print("""<text x="{x}" y="{y}" transform="rotate(90 {x},{y})" font-family="Courier" text-anchor="start" font-size="{fontsize}" fill="black" >{0:7.3f}</text> """.format(float(_)/float(factor), x = x, y = y + 30, fontsize = 12), file = io) # Drawing data points .. labels first for datasetnumber,dataset in enumerate(plot): r,g,b = self.colors[plotNumber][datasetnumber] style = self.styles[plotNumber][datasetnumber] if style[:5] == 'label': for pos,(mz,i) in enumerate(dataset): if resolved_mzRange[0] <= mz <= resolved_mzRange[1]: x = int(round(PADDING[3] + (mz - resolved_mzRange[0] ) * pixelpermz )) y = baseline y2 = baseline - (resolved_i_max * pixelper_i) x3 = round(x,-1) if style.upper().endswith('X1'): y3 = round(y+10) sign = 1 stroke_width = 0.50 else: y3 = round(y2) sign = -1 stroke_width = 0.40 layers = 15 while (x3,y3) in self.labelPos[plotNumber]: y3 = y3 + sign * 15 layers -= 1 if layers <= 0: break print("""<line x1="{0}" y1="{1}" x2="{0}" y2="{2}" style="stroke-opacity: 0.6;stroke:rgb({r},{g},{b});stroke-width:{o};" />""".format(x, y, y2, o=stroke_width, r=r,g=g,b=b),file=io) self.labelPos[plotNumber].add((x3,y3)) self.labelPos[plotNumber].add((x3+10,y3 + sign * 15)) print("""<text x="{x}" y="{y}" transform="rotate({angle} {x},{y})" font-family="Courier" text-anchor="start" font-size="{fontsize}" style="fill:rgb({r},{g},{b});" >{0}</text> """.format(i, x = x, y = y3, fontsize = 10, angle = sign * 45 , r = r, g = g, b = b), file = io) elif style[:5] == 'area': if self.normalizations[plotNumber] == True: maxI = max([ i for mz,i in dataset ]) else: maxI = None first = True for pos,(mz,i) in enumerate(sorted(dataset)): if resolved_mzRange[0] <= mz <= resolved_mzRange[1]: x = int(round(PADDING[3] + (mz - resolved_mzRange[0] ) * pixelpermz )) if self.normalizations[plotNumber] == True: y = baseline - ((float(i)/float(maxI)) * pixelper_i) else: y = baseline - (i * pixelper_i) if first: print('<path d="M{0} {2} L{0} {1}'.format(x,y,baseline), end=" ", file = io) first = False else: print('L{0} {1}'.format(x,y ), end=" ", file = io) if first == False: # aka we have seen some points ;) print(' L{0} {1} Z" style="fill:rgb({r},{g},{b}); fill-opacity:0.2; stroke:rgb({r},{g},{b}); stroke-width:1" />'.format(x,baseline,r=r,g=g,b=b), file = io) elif style[:8] == 'triangle': if self.normalizations[plotNumber] == True: maxI = max([ i for mz,i in dataset ]) else: maxI = None #first = True for pos,(mz,i) in enumerate(sorted(dataset)): if resolved_mzRange[0] <= mz <= resolved_mzRange[1]: x = int(round(PADDING[3] + (mz - resolved_mzRange[0] ) * pixelpermz )) xl = int(round(PADDING[3] + (mz-0.075 - resolved_mzRange[0] ) * pixelpermz )) xr = int(round(PADDING[3] + (mz+0.075 - resolved_mzRange[0] ) * pixelpermz )) if self.normalizations[plotNumber] == True: y = baseline - ((float(i)/float(maxI)) * pixelper_i) else: y = baseline - (i * pixelper_i) print('<path d="M{0} {4} L{2} {3} L{1} {4} Z" style="fill:rgb({r},{g},{b}); fill-opacity:0.4"/>'.format(xl,xr,x,y,baseline,r=r,g=g,b=b), file = io) #print('<line x1="{0}" y1="{1}" x2="{0}" y2="{2}" style="stroke:rgb({r},{g},{b});stroke-width:1"/>'.format(x,y,baseline,r=r,g=g,b=b), file = io) # Drawing data ... for datasetnumber,dataset in enumerate(plot): r,g,b = self.colors[plotNumber][datasetnumber] style = self.styles[plotNumber][datasetnumber] if style[:5] not in ['label', 'area', 'trian']: if len(dataset) == 0: continue if self.normalizations[plotNumber] == True: maxI = max([ i for mz,i in dataset ]) else: maxI = None for pos,(mz,i) in enumerate(dataset): if resolved_mzRange[0] <= mz <= resolved_mzRange[1]: x = int(round(PADDING[3] + (mz - resolved_mzRange[0] ) * pixelpermz )) if self.normalizations[plotNumber] == True: y = baseline - ((float(i)/float(maxI)) * pixelper_i) else: y = baseline - (i * pixelper_i) if style == 'circles': print("""<circle cx="{x}" cy="{y}" r="2" style="fill:rgb({r},{g},{b});" />""".format( x = x, y = y , r = r, g = g, b = b),file=io) elif style == 'emptycircles': print("""<circle cx="{x}" cy="{y}" r="4" style="stroke:rgb({r},{g},{b}); stroke-width:2; fill:none" />""".format( x = x, y = y , r = r, g = g, b = b),file=io) elif style == 'sticks': print("""<line x1="{0}" y1="{1}" x2="{0}" y2="{2}" style="stroke:rgb({r},{g},{b});stroke-width:0.5"/>""".format(x,y,baseline,r=r,g=g,b=b),file=io) elif style == 'squares': # this is experimental, the squares are not centred print("""<rect x="{0}" y="{1}" width="4" height="4" style="fill:rgb({r},{g},{b}); fill-opacity:0.8;" /> """.format(x-2,y-2,r=r,g=g,b=b),file=io) else: print("Style {0} not support yet ...".format(style),file=sys.stderr) exit(1) print("</svg>", file = io) print(""" </body> </html> """, file = io) return
if __name__ == '__main__': print(__doc__)