https://github.com/kpalin/EEL
Tip revision: c0cd936bedca30564088c8aeb8f43ebe4e05421f authored by Kimmo Palin on 07 November 2012, 13:23:26 UTC
Fixed a variable format clash which sometimes caused core dumps.
Fixed a variable format clash which sometimes caused core dumps.
Tip revision: c0cd936
pfm2ps.py
import sys
from cStringIO import StringIO
#
# $Log$
# Revision 1.1 2004/11/12 07:56:36 kpalin
# Initial check-in
#
#
# For making even sized sequence "logos" which are not sequence logos
# Legaleze for the postscript
##WebLogo License
##Copyright (c) 2002 Gavin E. Crooks, Gary Hon,
##John-Marc Chandonia and Steven E. Brenner
##Permission is hereby granted, free of charge, to
##any person obtaining a copy of this software and
##associated documentation files (the "Software"),
##to deal in the Software without restriction,
##including without limitation the rights to use,
##copy, modify, merge, publish, distribute,
##sublicense, and/or sell copies of the Software,
##and to permit persons to whom the Software is
##furnished to do so, subject to the following
##conditions:
##The above copyright notice and this permission
##notice shall be included in all copies or
##substantial portions of the Software.
##THE SOFTWARE IS PROVIDED "AS IS", WITHOUT
##WARRANTY OF ANY KIND, EXPRESS OR IMPLIED,
##INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
##MERCHANTABILITY, FITNESS FOR A PARTICULAR
##PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
##THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR
##ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER
##IN AN ACTION OF CONTRACT, TORT OR OTHERWISE,
##ARISING FROM, OUT OF OR IN CONNECTION WITH THE
##SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
##SOFTWARE.
##(This is the MIT Open Source License,
##http://www.opensource.org/licenses/mit-license.html)
class PfmFormat:
psProlog1="""%%!PS-Adobe-3.0 EPSF-3.0
%%%%Title: Sequence Logo:
%%%%Creator: pfm2ps
%%%%CreationDate:
%%%%BoundingBox: 0 0 %d 100
%%%%Pages: 0
%%%%DocumentFonts:
%%%%EndComments
"""
psProlog2="""
%% ---- CONSTANTS ----
/cmfactor 72 2.54 div def %% defines points -> cm conversion
/cm {cmfactor mul} bind def %% defines centimeters
%% ---- VARIABLES ----
/black [0 0 0] def
/red [0.8 0 0] def
/green [0 0.8 0] def
/blue [0 0 0.8] def
/yellow [1 0.71 .0] def
/purple [0.8 0 0.8] def
/orange [1 0.7 0] def
/charsPerLine %d def
/logoWidth 18 cm def
/logoHeight 5 cm def
/logoTitle () def
/yaxis false def
/yaxisLabel (bits) def
/yaxisBits 2 def %% bits
/yaxisTicBits 1 def
/xaxis false def
/xaxisLabel () def
/showEnds (d) def %% d: DNA, p: PROTEIN, -: none
/showFineprint false def
/fineprint (U. of Helsinki) def
/logoLines 1 def
/showingBox (n) def %%n s f
/shrinking false def
/shrink 0.5 def
/outline false def
/IbeamFraction 1 def
/IbeamGray 0.50 def
/IbeamLineWidth 0.5 def
/fontsize 12 def
/titleFontsize 14 def
/smallFontsize 6 def
/defaultColor black def
/charWidth %d def
"""
psProlog3="""
% Standard DNA/RNA color scheme
/colorDict <<
(G) orange
(T) blue
(C) green
(A) red
(U) black
>> def
% Standard DNA/RNA color scheme
% /colorDict <<
% (G) orange
% (T) red
% (C) blue
% (A) green
% (U) red
% >> def
% Standard Amino Acid colors
%/colorDict <<
% (G) green
% (S) green
% (T) green
% (Y) green
% (C) green
% (N) purple
% (Q) purple
% (K) blue
% (R) blue
% (H) blue
% (D) red
% (E) red
% (P) black
% (A) black
% (W) black
% (F) black
% (L) black
% (I) black
% (M) black
% (V) black
%>> def
% ---- DERIVED PARAMETERS ----
/leftMargin
fontsize 1.5 mul
def
/bottomMargin
fontsize 0.75 mul
% Add extra room for axis
xaxis {fontsize 1.75 mul add } if
xaxisLabel () eq {} {fontsize 0.75 mul add} ifelse
def
/topMargin
logoTitle () eq { 10 }{titleFontsize 4 add} ifelse
def
/rightMargin
%Add extra room if showing ends
showEnds (-) eq { fontsize}{fontsize 1.5 mul} ifelse
def
/yaxisHeight
logoHeight
bottomMargin sub
topMargin sub
def
/ticWidth fontsize 2 div def
/pointsPerBit yaxisHeight yaxisBits div def
/isBoxed
showingBox (s) eq
showingBox (f) eq or {
true
} {
false
} ifelse
def
/stackMargin 1 def
% Do not add space aroung characters if characters are boxed
/charRightMargin
isBoxed { 0.0 } {stackMargin} ifelse
def
/charTopMargin
isBoxed { 0.0 } {stackMargin} ifelse
def
%/charWidth
% logoWidth
% leftMargin sub
% rightMargin sub
% charsPerLine div
% charRightMargin sub
%def
/charWidth4 charWidth 4 div def
/charWidth2 charWidth 2 div def
/stackWidth
charWidth charRightMargin add
def
/numberFontsize
fontsize charWidth lt {fontsize}{charWidth} ifelse
def
% movements to place 5'/N and 3'/C symbols
/leftEndDeltaX fontsize neg def
/leftEndDeltaY fontsize 1.5 mul neg def
/rightEndDeltaX fontsize 0.25 mul def
/rightEndDeltaY leftEndDeltaY def
% Outline width is proporional to charWidth,
% but no less that 1 point
/outlinewidth
charWidth 32 div dup 1 gt {}{pop 1} ifelse
def
% ---- PROCEDURES ----
/StartLogo {
% Save state
save
gsave
% Print Logo Title, top center
gsave
SetTitleFont
logoWidth 2 div
logoTitle
stringwidth pop 2 div sub
logoHeight logoLines mul
titleFontsize sub
moveto
logoTitle
show
grestore
% Print X-axis label, bottom center
gsave
SetStringFont
logoWidth 2 div
xaxisLabel stringwidth pop 2 div sub
fontsize 3 div
moveto
xaxisLabel
show
grestore
% Show Fine Print
showFineprint {
gsave
SetSmallFont
logoWidth
fineprint stringwidth pop sub
smallFontsize sub
smallFontsize 3 div
moveto
fineprint show
grestore
} if
% Move to lower left corner of last line, first stack
leftMargin bottomMargin translate
% Move above first line ready for StartLine
0 logoLines logoHeight mul translate
SetLogoFont
} bind def
/EndLogo {
grestore
showpage
restore
} bind def
/StartLine{
% move down to the bottom of the line:
0 logoHeight neg translate
gsave
yaxis { MakeYaxis } if
xaxis { ShowLeftEnd } if
} bind def
/EndLine{
xaxis { ShowRightEnd } if
grestore
} bind def
/MakeYaxis {
gsave
stackMargin neg 0 translate
ShowYaxisBar
ShowYaxisLabel
grestore
} bind def
/ShowYaxisBar {
gsave
SetStringFont
/str 10 string def % string to hold number
/smallgap stackMargin 2 div def
% Draw first tic and bar
gsave
ticWidth neg 0 moveto
ticWidth 0 rlineto
0 yaxisHeight rlineto
stroke
grestore
% Draw the tics
% initial increment limit proc for
0 yaxisTicBits yaxisBits abs %cvi
{/loopnumber exch def
% convert the number coming from the loop to a string
% and find its width
loopnumber 10 str cvrs
/stringnumber exch def % string representing the number
stringnumber stringwidth pop
/numberwidth exch def % width of number to show
/halfnumberheight
stringnumber CharBoxHeight 2 div
def
numberwidth % move back width of number
neg loopnumber pointsPerBit mul % shift on y axis
halfnumberheight sub % down half the digit
moveto % move back the width of the string
ticWidth neg smallgap sub % Move back a bit more
0 rmoveto % move back the width of the tic
stringnumber show
smallgap 0 rmoveto % Make a small gap
% now show the tic mark
0 halfnumberheight rmoveto % shift up again
ticWidth 0 rlineto
stroke
} for
grestore
} bind def
/ShowYaxisLabel {
gsave
SetStringFont
% How far we move left depends on the size of
% the tic labels.
/str 10 string def % string to hold number
yaxisBits yaxisTicBits div cvi yaxisTicBits mul
str cvs stringwidth pop
ticWidth 1.5 mul add neg
yaxisHeight
yaxisLabel stringwidth pop
sub 2 div
translate
90 rotate
0 0 moveto
yaxisLabel show
grestore
} bind def
/StartStack { % <stackNumber> startstack
xaxis {MakeNumber}{pop} ifelse
gsave
} bind def
/EndStack {
grestore
stackWidth 0 translate
} bind def
% Draw a character whose height is proportional to symbol bits
/MakeSymbol{ % charbits character MakeSymbol
gsave
/char exch def
/bits exch def
/bitsHeight
bits pointsPerBit mul
def
/charHeight
bitsHeight charTopMargin sub
dup
0.0 gt {}{pop 0.0} ifelse % if neg replace with zero
def
charHeight 0.0 gt {
char SetColor
charWidth charHeight char ShowChar
showingBox (s) eq { % Unfilled box
0 0 charWidth charHeight false ShowBox
} if
showingBox (f) eq { % Filled box
0 0 charWidth charHeight true ShowBox
} if
} if
grestore
0 bitsHeight translate
} bind def
/ShowChar { % <width> <height> <char> ShowChar
gsave
/tc exch def % The character
/ysize exch def % the y size of the character
/xsize exch def % the x size of the character
/xmulfactor 1 def
/ymulfactor 1 def
% if ysize is negative, make everything upside down!
ysize 0 lt {
% put ysize normal in this orientation
/ysize ysize abs def
xsize ysize translate
180 rotate
} if
shrinking {
xsize 1 shrink sub 2 div mul
ysize 1 shrink sub 2 div mul translate
shrink shrink scale
} if
% Calculate the font scaling factors
% Loop twice to catch small correction due to first scaling
2 {
gsave
xmulfactor ymulfactor scale
ysize % desired size of character in points
tc CharBoxHeight
dup 0.0 ne {
div % factor by which to scale up the character
/ymulfactor exch def
} % end if
{pop pop}
ifelse
xsize % desired size of character in points
tc CharBoxWidth
dup 0.0 ne {
div % factor by which to scale up the character
/xmulfactor exch def
} % end if
{pop pop}
ifelse
grestore
} repeat
% Adjust horizontal position if the symbol is an I
tc (I) eq {
charWidth 2 div % half of requested character width
tc CharBoxWidth 2 div % half of the actual character
sub 0 translate
% Avoid x scaling for I
/xmulfactor 1 def
} if
% ---- Finally, draw the character
newpath
xmulfactor ymulfactor scale
% Move lower left corner of character to start point
tc CharBox pop pop % llx lly : Lower left corner
exch neg exch neg
moveto
outline { % outline characters:
outlinewidth setlinewidth
tc true charpath
gsave 1 setgray fill grestore
clip stroke
} { % regular characters
tc show
} ifelse
grestore
} bind def
/ShowBox { % x1 y1 x2 y2 filled ShowBox
gsave
/filled exch def
/y2 exch def
/x2 exch def
/y1 exch def
/x1 exch def
newpath
x1 y1 moveto
x2 y1 lineto
x2 y2 lineto
x1 y2 lineto
closepath
clip
filled {
fill
}{
0 setgray stroke
} ifelse
grestore
} bind def
/MakeNumber { % number MakeNumber
gsave
SetNumberFont
stackWidth 0 translate
90 rotate % rotate so the number fits
dup stringwidth pop % find the length of the number
neg % prepare for move
stackMargin sub % Move back a bit
charWidth (0) CharBoxHeight % height of numbers
sub 2 div %
moveto % move back to provide space
show
grestore
} bind def
/Ibeam{ % heightInBits Ibeam
gsave
% Make an Ibeam of twice the given height in bits
/height exch pointsPerBit mul def
/heightDRAW height IbeamFraction mul def
IbeamLineWidth setlinewidth
IbeamGray setgray
charWidth2 height neg translate
ShowIbar
newpath
0 0 moveto
0 heightDRAW rlineto
stroke
newpath
0 height moveto
0 height rmoveto
currentpoint translate
ShowIbar
newpath
0 0 moveto
0 heightDRAW neg rlineto
currentpoint translate
stroke
grestore
} bind def
/ShowIbar { % make a horizontal bar
gsave
newpath
charWidth4 neg 0 moveto
charWidth4 0 lineto
stroke
grestore
} bind def
/ShowLeftEnd {
gsave
SetStringFont
leftEndDeltaX leftEndDeltaY moveto
showEnds (d) eq {(5) show ShowPrime} if
showEnds (p) eq {(N) show} if
grestore
} bind def
/ShowRightEnd {
gsave
SetStringFont
rightEndDeltaX rightEndDeltaY moveto
showEnds (d) eq {(3) show ShowPrime} if
showEnds (p) eq {(C) show} if
grestore
} bind def
/ShowPrime {
gsave
SetPrimeFont
(\242) show
grestore
} bind def
/SetColor{ % <char> SetColor
dup colorDict exch known {
colorDict exch get aload pop setrgbcolor
} {
pop
defaultColor aload pop setrgbcolor
} ifelse
} bind def
% define fonts
/SetTitleFont {/Times-Bold findfont titleFontsize scalefont setfont} bind def
%/SetLogoFont {/Helvetica-Narrow-Bold findfont charWidth scalefont setfont} bind def
/SetLogoFont {/Times-Bold findfont charWidth scalefont setfont} bind def
/SetStringFont{/Helvetica-Bold findfont fontsize scalefont setfont} bind def
/SetPrimeFont {/Symbol findfont fontsize scalefont setfont} bind def
/SetSmallFont {/Helvetica findfont smallFontsize scalefont setfont} bind def
/SetNumberFont {
/Helvetica-Bold findfont
numberFontsize
scalefont
setfont
} bind def
%Take a single character and return the bounding box
/CharBox { % <char> CharBox <lx> <ly> <ux> <uy>
gsave
newpath
0 0 moveto
% take the character off the stack and use it here:
true charpath
flattenpath
pathbbox % compute bounding box of 1 pt. char => lx ly ux uy
% the path is here, but toss it away ...
grestore
} bind def
% The height of a characters bounding box
/CharBoxHeight { % <char> CharBoxHeight <num>
CharBox
exch pop sub neg exch pop
} bind def
% The width of a characters bounding box
/CharBoxWidth { % <char> CharBoxHeight <num>
CharBox
pop exch pop sub neg
} bind def
% Deprecated names
/startstack {StartStack} bind def
/endstack {EndStack} bind def
/makenumber {MakeNumber} bind def
/numchar { MakeSymbol } bind def
%%%%EndProlog
%%%%Page: 1 1
StartLogo
StartLine %% line number 1
"""
psEpilog="""
EndLine
EndLogo
%%EOF
"""
def __init__(self,fname=None,sortOrder=1):
self.fname=fname
self.charwidth=40
if self.fname:
self.pfm=self.readpfm(self.fname)
# From top
# 0-A,C,G,T
# 1-Larger on Top
# 2-Smaller on Top
self.sortOrder=sortOrder
def readpfm(self,fname):
return [map(int,x.split()) for x in open(fname).readlines()]
def makeStacks(self,pfm=None):
if pfm==None:
pfm=self.pfm
self.stacks=len(pfm[0])
stackCodes=["A","C","G","T"]
self.stackStrs=[]
for i in range(self.stacks):
tot=0
for j in range(len(pfm)):
tot=tot+pfm[j][i]
sStr=[]
for j in range(len(pfm)):
if pfm[j][i]>0:
height=pfm[j][i]*1.0/tot
char=stackCodes[j]
sStr.append((height,char," %f (%s) numchar\n"%(height,char)))
self.stackStrs.append(sStr)
return self.stackStrs
def setSort(self,sortOrder):
self.sortOrder=sortOrder
def sort(self):
for i in range(self.stacks):
if self.sortOrder==0:
self.stackStrs[i].sort(lambda x,y:cmp(y[1],x[1]))
elif self.sortOrder==1:
self.stackStrs[i].sort(lambda x,y:cmp((x[0],-ord(x[1])),(y[0],-ord(y[1]))))
elif self.sortOrder==2:
self.stackStrs[i].sort(lambda x,y:cmp((y[0],ord(y[1])),(x[0],ord(x[1]))))
def stacksStr(self,s=StringIO()):
self.sort()
for (i,stack) in zip(range(self.stacks),self.stackStrs):
s.write("(%d) startstack\n"%(i+1))
for height,char,psString in stack:
s.write(psString)
s.write("endstack\n\n")
return s
def __str__(self):
s=StringIO()
if not hasattr(self,"stacks"):
self.makeStacks()
s.write(self.psProlog1%((self.stacks+1)*self.charwidth))
s.write(self.psProlog2%(self.stacks,self.charwidth))
s.write(self.psProlog3)
self.stacksStr(s)
s.write(self.psEpilog)
return s.getvalue()
if __name__=="__main__":
if len(sys.argv)<2:
print "usage: python pfm2ps.py TF/joku.pfm [sortcode]"
sys.exit(1)
sortOrder=1
if len(sys.argv)>=3:
sortOrder=int(sys.argv[2])
pfm=PfmFormat(sys.argv[1],sortOrder)
print str(pfm)
#print "\n".join(makeStacks(pfm))