summaryrefslogtreecommitdiff
path: root/macros/latex/contrib/cpssp/cpssp
diff options
context:
space:
mode:
Diffstat (limited to 'macros/latex/contrib/cpssp/cpssp')
-rwxr-xr-xmacros/latex/contrib/cpssp/cpssp412
1 files changed, 0 insertions, 412 deletions
diff --git a/macros/latex/contrib/cpssp/cpssp b/macros/latex/contrib/cpssp/cpssp
deleted file mode 100755
index 40d825f05a..0000000000
--- a/macros/latex/contrib/cpssp/cpssp
+++ /dev/null
@@ -1,412 +0,0 @@
-#!/usr/bin/python
-# CPSSP -- Compare Protein Secondary Structure Predictions
-# v1.0 20090606
-# Copyright (C) 2009 by Wolfgang Skala
-#
-# This work may be distributed and/or modified under the
-# conditions of the LaTeX Project Public License, either version 1.3
-# of this license or (at your option) any later version.
-# The latest version of this license is in
-# http://www.latex-project.org/lppl.txt
-# and version 1.3 or later is part of all distributions of LaTeX
-# version 2005/12/01 or later.
-
-import getopt, sys
-
-
-
-### 1. FUNCTIONS
-
-def readFasta(filename, alphabet):
- # reads a FASTA file containing one or several sequences/alignments/predictions
- # input: * (filename) the name of the FASTA file
- # * (alphabet) the allowed characters in the file, e.g. amino acids
- # output: * (names) a list of strings corresponding to the sequence names
- # * (allseqs) a list containing the sequences;
- # each sequence is a list of characters, where each character represents a
- # residue, gap or secondary structure element
- allseqs = []
- names = []
- seq = []
- f = file(filename, "r")
- for line in f.readlines():
- if line[0] == ">":
- if seq != []:
- allseqs.append(seq)
- seq = []
- names.append(line[1:-1])
- else:
- for c in line:
- c = c.upper()
- if c in alphabet:
- seq.append(c)
- if seq != []:
- allseqs.append(seq)
- f.close()
- return names, allseqs
-
-
-def removeGaps(allSeqs):
- # removes gaps from a list of sequences
- # if a position in each sequence is occupied exclusively by gaps
- # input: (allSeqs) a list of sequences as returned by readFasta()
- # output: same format as input with gaps removes
- gapPositions = []
- for i in range(len(allSeqs[0])-1, -1, -1):
- onlyGaps = True
- for j in range(len(allSeqs)):
- if allSeqs[j][i] not in GAP_SYMBOLS:
- onlyGaps = False
- if onlyGaps:
- for j in range(len(allSeqs)):
- allSeqs[j].pop(i)
- return allSeqs
-
-
-def commonSyntax(allStrucs):
- # make all structures a common syntax, i.e. a coil is represented by "C" and not by "-" or a space
- # input: (allStrucs) list of structures as returned by readFasta
- # output: same format as input
- for i in range(len(allStrucs)):
- for j in range(len(allStrucs[i])):
- if allStrucs[i][j] == "-" or allStrucs[i][j] == " ":
- allStrucs[i][j] = "C"
- return allStrucs
-
-
-def addGaps(allSeqs, allStrucs):
- # add gaps to the secondary structures so that they correspond to the gapped sequences
- # input: * (allSeqs) list of (degapped) sequences as returned by removeGaps() OR None,
- # indicating that no gaps should be added
- # * (allStrucs) list of structures as returned by readFasta/commonSyntax
- # output: list of gapped structures (same format as allStrucs)
- resultStrucs = []
- if allSeqs == None:
- for struc in allStrucs:
- curStruc = ""
- for res in struc:
- curStruc += res
- resultStrucs.append(curStruc)
- else:
- i = 0
- for seq in allSeqs:
- k = 0
- curStruc = ""
- for j in range(len(seq)):
- if seq[j] not in GAP_SYMBOLS:
- curStruc += allStrucs[i][k]
- k += 1
- else:
- curStruc += "-"
- resultStrucs.append(curStruc)
- i += 1
- return resultStrucs
-
-
-def breakLines(allStrucs, n):
- # break the structures into lines according to the number of residues per line specified by
- # the user if a line ends with a sheet (E) and the next line starts with a sheet, change the end
- # letter to "e" which indicates that no arrowhead should be drawn in the graphical representation
- # input: (allStrucs) list of gapped structures as returned by addGaps()
- # (n) residues per line
- # output: list of structures; each structure is list of 3-tuplets containing (1) a string which
- # specifies the residues on the line, (2) the number of the first and (3) the number of
- # the last residue on the line
- resultStrucs = []
- for struc in allStrucs:
- curStruc = []
- for i in range(len(struc) / n + 1):
- if i*n != len(struc):
- s = struc[i*n:(i+1)*n]
- if s[-1] == "E" and (i+1)*n < len(struc):
- if struc[(i+1)*n] == "E":
- curStruc.append(s[:-1] + "e")
- else:
- curStruc.append(s)
- else:
- curStruc.append(s)
- resultStrucs.append(curStruc)
- for struc in resultStrucs:
- startRes = 0
- endRes = 0
- for i in range(len(struc)):
- for j in range(len(struc[i])):
- if struc[i][j] != "-":
- endRes += 1
- struc[i] = (struc[i], startRes, endRes)
- startRes = endRes
- return resultStrucs
-
-
-def makeTikzDraw(ssType, block, line, start, end):
- # compose a TikZ command which draws a secondary structure element
- # input: * (ssType) secondary structure type (C, coil; H, helix; E, sheet; e, sheet at end of
- # the line; -, gap)
- # * (block) the current sequence block
- # * (line) the current line
- # * (start) the start position
- # * (end) the end position
- # output: a string containing the TikZ command
- result = "\t\t"
- if ssType == "B":
- result += "\\cpsspBridge{-"
- elif ssType == "C":
- result += "\\cpsspCoil{-"
- elif ssType == "E":
- result += "\\cpsspSheet{-"
- elif ssType == "e":
- result += "\\cpsspSheetT{-"
- elif ssType == "G":
- result += "\\cpsspThreeTenHelix{-"
- elif ssType == "H":
- result += "\\cpsspAlphaHelix{-"
- elif ssType == "I":
- result += "\\cpsspPiHelix{-"
- elif ssType == "S":
- result += "\\cpsspBend{-"
- elif ssType == "T":
- result += "\\cpsspTurn{-"
- else:
- result += "\\cpsspGap{-"
-
- result += str(block * (blockDistance + nStruc * lineDistance) + line * lineDistance) + "}{"
- result += str(lineIndent + resWidth * start) + "}{"
- result += str(lineIndent + resWidth * end) + "}\n"
- return result
-
-
-def makeTikzLabel(text, block, line):
- # compose a TikZ command which draws a label
- # input: * (text) the label text
- # * (block) the current sequence block
- # * (line) the current line
- # output: a string containing the TikZ command
- result = "\t\\cpsspLabel{-"
- result += str(block * (blockDistance + nStruc * lineDistance) + line * lineDistance) + "}{"
- result += text + "}\n"
- return result
-
-
-def makeTikzRes(number, block, line, isStart, pos=None):
- # compose a TikZ command which draws the number of the first residue in the line
- # input: * (number) the residue number
- # * (block) the current sequence block
- # * (line) the current line
- # * (isStart) True if the start residue number is to be drawn, False otherwise
- # * (pos) unused for the start residue; for the end residue, it indicates the x position
- result = "\t\t"
- if isStart:
- result +="\\cpsspStartRes{-"
- result += str(block * (blockDistance + nStruc * lineDistance) + line * lineDistance) + "}{"
- result += str(lineIndent) + "}{"
- else:
- result +="\\cpsspEndRes{-"
- result += str(block * (blockDistance + nStruc * lineDistance) + line * lineDistance) + "}{"
- result += str(lineIndent + resWidth * pos) + "}{"
- result += str(number) + "}\n"
- return result
-
-
-def usage():
- # print usage of the program
- print """CPSSP -- Compare Protein Secondary Structure Prediction v1.0
-Usage: cpssp
--h or --help prints this message
--s or --sequence-file (FASTA file containing the sequences)
--u or --structure-file (FASTA file containing the structures; mandatory)
--o or --output-file (outout filename without extension and numbering)
--w or --image-width (width of the image in cm)
--t or --image-height (maximal height of an image in cm)
--i or --line-indent (indentation at the beginning of the line in cm)
--r or --residues-per-line (number of residues per line)
--l or --line-distance (distance between lines in cm)
--b or --block-distance (distance between blocks in cm)"""
-
-
-def version():
- # print the program version
- print """CPSSP 1.0
-Copyright (C) 2009 Wolfgang Skala
-License LPPL v1.3c: The LaTeX project public license version 1.3c <http://www.latex-project.org/lppl.txt>
-This is free software: you are free to change and redistribute it.
-There is NO WARRANTY, to the extent permitted by law."""
-
-
-### 2. CONSTANTS AND VARIABLES
-
-AMINO_ACIDS = 'ARNDCQEGHILKMFPSTWYV-.' # the characters allowed in the FASTA file (AAs, gaps, structures)
-SS_ELEMENTS = "BCEGHIST- "
-GAP_SYMBOLS = "-." # possible gap symbols
-
-sequenceFile = None
-structureFile = None
-imageWidth = 15 # total line width (in cm)
-imageHeight = 20 # maximal height of the image (in cm); it will be split into separate files
- # if its natural height exceeds this value; 0 indicates an arbitrary height
-lineIndent = 2.5 # indentation at the left pof each line (in cm)
-resPerLine = 50 # number of residues per line
-lineDistance = .5 # distance between sequences within one line (in cm)
-blockDistance = 1 # distance between sequence blocks (in cm)
-outputFile = "cpsspresult"
-
-
-
-
-### 3. MAIN PART
-
-# process command line options
-try:
- opts, args = getopt.getopt(sys.argv[1:],
- "vhs:u:w:t:r:i:l:b:o:",
- ["version", "help", "sequence-file=", "structure-file=", "image-width=", "image-height=",
- "residues-per-line=", "line-indent=", "line-distance=", "block-distance=", "output-file="])
-except getopt.GetoptError:
- usage()
- sys.exit(2)
-
-for opt, arg in opts:
- if opt in ("-h", "--help"):
- usage()
- sys.exit()
- elif opt in ("-v", "--version"):
- version()
- sys.exit()
- elif opt in ("-s", "--sequence-file"):
- sequenceFile = arg
- elif opt in ("-u", "--structure-file"):
- structureFile = arg
- elif opt in ("-o", "--output-file"):
- outputFile = arg
- elif opt in ("-w", "--image-width"):
- try:
- imageWidth = float(arg)
- except ValueError:
- print "Invalid image width."
- sys.exit(1)
- elif opt in ("-t", "--image-height"):
- try:
- imageHeight = float(arg)
- except ValueError:
- print "Invalid image height."
- sys.exit(1)
- elif opt in ("r", "--residues-per-line"):
- try:
- resPerLine = int(arg)
- except ValueError:
- print "Invalid number of residues per line."
- sys.exit(1)
- elif opt in ("-i", "--line-indent"):
- try:
- lineIndent = float(arg)
- except ValueError:
- print "Invalid line indentation."
- sys.exit(1)
- elif opt in ("-l", "--line-distance"):
- try:
- lineDistance = float(arg)
- except ValueError:
- print "Invalid line distance."
- sys.exit(1)
- elif opt in ("-b", "--block-distance"):
- try:
- blockDistance = float(arg)
- except ValueError:
- print "Invalid block distance."
- sys.exit(1)
-
-if structureFile == None:
- usage()
- sys.exit(2)
-elif sequenceFile == None:
- # compare the predictions from multiple programs for a single protein
- # open FASTA file
- try:
- seqNames, structures = readFasta(structureFile, SS_ELEMENTS)
- except IOError as error:
- print "Could not open '" + error.filename + "'."
- sys.exit(1)
-
- # process structures
- try:
- structures = commonSyntax(structures)
- structures = addGaps(None, structures)
- brokenStructures = breakLines(structures, resPerLine)
- except IndexError:
- print "The structures seem to differ in length."
- sys.exit(1)
-else:
- # compare the predictions from a single program for multiple proteins
- # open FASTA files
- try:
- seqNames, sequences = readFasta(sequenceFile, AMINO_ACIDS)
- strucNames, structures = readFasta(structureFile, SS_ELEMENTS)
- except IOError as error:
- print "Could not open '" + error.filename + "'."
- sys.exit(1)
-
- # process sequences and structures read from the files
- try:
- sequences = removeGaps(sequences)
- structures = commonSyntax(structures)
- structures = addGaps(sequences, structures)
- brokenStructures = breakLines(structures, resPerLine)
- except IndexError:
- print "The sequences and structures seem to differ in length."
- sys.exit(1)
-
-
-# now for the common part
-# calculate or initiate some variables
-resWidth = float(imageWidth - lineIndent) / resPerLine # width of a single residue (in cm)
-nStruc = len(structures) # number of structures
-nBlocks = len(brokenStructures[0]) # number of blocks
-if imageHeight == 0:
- blocksPerImage = nBlocks # blocks per image (output file)
-else:
- blocksPerImage = int((imageHeight + blockDistance) / (nStruc * lineDistance + blockDistance))
-
-tikzCommands = [] # list of strings where each string contains all TikZ commands for an image
-for i in range(nBlocks / blocksPerImage + (1 if nBlocks % blocksPerImage != 0 else 0)):
- tikzCommands.append("")
-
-# determine the appropriate commands
-curLine = 0
-for struc in brokenStructures:
- curBlock = 0
- curImage = 0
- for line, startRes, endRes in struc:
- tikzCommands[curImage] += makeTikzLabel(seqNames[curLine], curBlock, curLine)
- curType = ""
- for i in range(len(line)):
- if i == 0:
- curType = line[i]
- startPos = i
- if line[i] != "-":
- tikzCommands[curImage] += makeTikzRes(startRes + 1, curBlock, curLine, True)
- else:
- tikzCommands[curImage] += makeTikzRes(startRes, curBlock, curLine, True)
- if i < len(line) - 1:
- if line[i+1].upper() != line[i]:
- tikzCommands[curImage] += makeTikzDraw(curType, curBlock, curLine, startPos, i + 1)
- curType = line[i+1]
- startPos = i + 1
- else:
- if line[i] == "e":
- tikzCommands[curImage] += makeTikzDraw("e", curBlock, curLine, startPos, i + 1)
- else:
- tikzCommands[curImage] += makeTikzDraw(curType, curBlock, curLine, startPos, i + 1)
- tikzCommands[curImage] += makeTikzRes(endRes, curBlock, curLine, False, i + 1)
- curBlock += 1
- if curBlock % blocksPerImage == 0:
- curImage += 1
- curBlock = 0
- curLine += 1
-
-# write the output files
-try:
- for i in range(len(tikzCommands)):
- f = file(outputFile + str(i) + ".tex", "w")
- f.write(tikzCommands[i])
- f.close()
-except IOError as error:
- print "Error while writing '" + error.filename + "'."