[med-svn] [Git][med-team/suitename][master] 8 commits: Code says its version 0.4.130509
Andreas Tille (@tille)
gitlab at salsa.debian.org
Fri Jan 28 16:25:11 GMT 2022
Andreas Tille pushed to branch master at Debian Med / suitename
Commits:
d2a32f6d by Andreas Tille at 2022-01-28T16:57:11+01:00
Code says its version 0.4.130509
- - - - -
1fbd9b78 by Andreas Tille at 2022-01-28T17:15:44+01:00
Upstream moved to Python and the original suitename code was moved - see README.Debian. This is the last commit with the old C code
- - - - -
372f2f07 by Andreas Tille at 2022-01-28T17:16:50+01:00
New upstream version 0.4.130509+git20210223.ebb1325
- - - - -
dfe960f2 by Andreas Tille at 2022-01-28T17:16:50+01:00
Update upstream source from tag 'upstream/0.4.130509+git20210223.ebb1325'
Update to upstream version '0.4.130509+git20210223.ebb1325'
with Debian dir 0890d72e8132d46e71547ab84327bd4c25334a96
- - - - -
9515f562 by Andreas Tille at 2022-01-28T17:20:24+01:00
Add README.Debian to explain status of this code
- - - - -
6c104fc8 by Andreas Tille at 2022-01-28T17:20:49+01:00
Update copyright
- - - - -
0c17255d by Andreas Tille at 2022-01-28T17:21:09+01:00
routine-update: Standards-Version: 4.6.0
- - - - -
e3d67f08 by Andreas Tille at 2022-01-28T17:21:52+01:00
routine-update: Ready to upload to unstable
- - - - -
13 changed files:
- + CMakeLists.txt
- + README.md
- + debian/README.Debian
- debian/changelog
- debian/control
- debian/copyright
- debian/watch
- + python/suitename.py
- + python/suitenamedefs.py
- + python/suiteninit.py
- + python/suiteninput.py
- + python/suitenout.py
- + python/suitenutil.py
Changes:
=====================================
CMakeLists.txt
=====================================
@@ -0,0 +1,29 @@
+cmake_minimum_required(VERSION 3.10.0)
+project(suitename)
+
+if (WIN32)
+ add_definitions(-D_CRT_SECURE_NO_WARNINGS)
+endif ()
+
+set(suitename_SOURCES
+ suitename.c
+ suiteninit.c
+ suitenscrt.c
+ suiteninpt.c
+ suitenout.c
+ suitenutil.c
+)
+
+set(suitename_HEADERS
+ suitename.h
+ suiteninit.h
+ suitenscrt.h
+ suiteninpt.h
+ suitenout.h
+ suitenutil.h
+)
+
+add_executable(suitename ${suitename_SOURCES} ${suitename_HEADERS})
+if (NOT WIN32)
+ target_link_libraries(suitename m)
+endif()
=====================================
README.md
=====================================
@@ -0,0 +1,21 @@
+# suitename
+Suitename - a tool for classifying backbone "suite" conformations (linkages between ribose sugars) in an RNA molecule.
+
+Suitename takes an input file as its first argument, or if none provided, reads the standard input. It writes its results to standard output. It classifies the conformation of each suite into one of several dozen predefined clusters that have been determined by years of study.
+
+Two forms of input are supported:
+ 1. A list of the 6 dihedral angles $/alpha/, $/beta/,$/gamma/,$/delta/, $/epsilon/,$/zeta/) for each residue of the RNA molecule. Suitename will re-parse the dihedral angles in the residues to obtain the 7 dihedral angles in each suite ($/delta/-1, $/epsilon/-1, $/zeta/-1, $/alpha/, $/beta/, $/gamma/, $/delta/), and then operate on those. This is the default input.
+
+ Each line of this format describes one residue, with fields separated by colons. The first several (default 6) are ID information, the remainder are angles. Sample:
+1: A: 7: : : U:-75.533:-154.742:48.162:80.895:-148.423:-159.688
+
+ 2. A kinemage file providing a list of 7 or 9 dihedral angles in each suite. Mark this by using the --suitein command line flag. 9 angles include $/chi/-1 and $/chi/.
+
+Three forms of output are supported:
+ 1. A text report, showing the classification of each suite into a cluster or outlier, and the neatness of fit ("suiteness") of the suite into that cluster. Suiteness is the cosine of the normalized distance of a suite datapoint from the power 3 hyperellipsoid boundary of the cluster in toward its center. A statistical summary is included at the end. This is the default output format.
+ 2. A kinemage file, which will display a 7D data point for each suite in the data. Points are color-coded according to the clusters to which they have been assigned. Each cluster is displayed as a colored ring surrounding its defined center. Specify this by using the --kinemage command line flag.
+ 3. A brief string showing only the cluster assignments. It consists of three characters per suite - base identity (uc) and 2-character number-letter name of the suite cluster (e.g., C1aG1gU1aA1aA1cG). Specify this by using the --string command line flag.
+
+(Hyperellipsoids: http://www.cs.utah.edu/dept/old/texinfo/glibc-manual-0.02/library_17.html)
+
+Many other command line options are available; type suitename --help to display them.
=====================================
debian/README.Debian
=====================================
@@ -0,0 +1,16 @@
+suitename for Debian
+====================
+
+Obsolete
+--------
+
+This package contains obsolete code that is no longer maintained as of August 2021.
+SuiteName development has been moved into the CCTBX code base and it can
+be run from a CCTBX installation as phenix.suitename or molprobity.suitename.
+Its location within that code base is here:
+ https://github.com/cctbx/cctbx_project/mmtbx/suitename).
+
+This package contains the last commit of the original C implementation of SuiteName
+along with build files to compile into a running executable.
+
+ -- Andreas Tille <tille at debian.org> Fri, 28 Jan 2022 16:56:19 +0100
=====================================
debian/changelog
=====================================
@@ -1,3 +1,11 @@
+suitename (0.4.130509+git20210223.ebb1325-1) unstable; urgency=medium
+
+ * Code says its version 0.4.130509
+ * Add README.Debian to explain status of this code
+ * Standards-Version: 4.6.0 (routine-update)
+
+ -- Andreas Tille <tille at debian.org> Fri, 28 Jan 2022 17:21:24 +0100
+
suitename (0.3.070919+git20180613.ebb1325-2) unstable; urgency=medium
* Team Upload.
=====================================
debian/control
=====================================
@@ -5,7 +5,7 @@ Uploaders: Andreas Tille <tille at debian.org>,
Section: science
Priority: optional
Build-Depends: debhelper-compat (= 13)
-Standards-Version: 4.5.0
+Standards-Version: 4.6.0
Vcs-Browser: https://salsa.debian.org/med-team/suitename
Vcs-Git: https://salsa.debian.org/med-team/suitename.git
Homepage: http://kinemage.biochem.duke.edu/software/suitename.php
=====================================
debian/copyright
=====================================
@@ -7,7 +7,7 @@ Copyright: 2007 David C. Richardson
License: suitename
Files: debian/*
-Copyright: 2015 Andreas Tille <tille at debian.org>
+Copyright: 2015-2022 Andreas Tille <tille at debian.org>
License: suitename
License: suitename
=====================================
debian/watch
=====================================
@@ -1,9 +1,9 @@
version=4
+opts=dversionmangle=s/.*/0.No-Release/ \
+ https://people.debian.org/~eriberto/ FakeWatchNoUpstreamReleaseForThisPackage-(\d\S+)\.gz
-opts="mode=git,pretty=0.3.070919+git%cd.%h" \
- https://github.com/rlabduke/suitename.git HEAD
+# Last version in C - see README.Debian
+#opts="mode=git,pretty=0.4.130509+git%cd.e37b832" \
+# https://github.com/rlabduke/suitename.git e37b832
-# http://kinemage.biochem.duke.edu/software/suitename.php .*/downloads/software/suitename/suitename.([\d.]+)\.src\.tgz
-
-# if there would be release tags ...
-# https://github.com/rlabduke/suitename/releases .*/archive/#PREFIX#(\d[\d.-]+)\.(?:tar(?:\.gz|\.bz2)?|tgz)
+# https://github.com/rlabduke/suitename.git HEAD ### HEAD is now Python code (see README.Debian
=====================================
python/suitename.py
=====================================
@@ -0,0 +1,485 @@
+import suiteninit, suitenout
+from suitenamedefs import Suite, Residue, Bin, Cluster, Issue, failMessages
+from suiteninit import args, bins, MAX_CLUSTERS
+from suiteninit import normalWidths, satelliteWidths
+from suiteninput import readResidues, readKinemageFile, buildSuites
+from suitenutil import hyperEllipsoidDistance
+
+import sys, os
+import numpy as np
+from math import cos, pi
+
+# suitename.py
+# ***************************************************************
+# NOTICE: This is free software and the source code is freely
+# available. You are free to redistribute or modify under the
+# conditions that (1) this notice is not removed or modified
+# in any way and (2) any modified versions of the program are
+# also available for free.
+# ** Absolutely no Warranty **
+# Copyright (C) 2007 David C. Richardson
+# ***************************************************************
+
+# 0.2.070524 preserve chi-1 and chi, so could preserve eta, theta
+# 0.3.070525 general read dangle record for, e.g., eta, theta
+# 0.3.070628 triage reports zeta-1, epsilon-1, delta-1,... Ltriage codes
+# 0.3.070803 notes: rearranged suitenhead.h/janesviews ...
+# put something in to say what veiws mean
+# what are masters e and d ????
+# 0.3.070919 3g wannabe (tRNA TpseudoUC loop)
+# 0.3.110606 range of delta updated by S.J.
+# 01/07/2014 S.J. updated so that it can take input with alternate conformations, *nd by default will calculate the suite for altA
+# 09/18/2014 S.J. updated so that suitename will ignore DNA residues
+
+
+version = "suitename.0.6.012521"
+dbCounter = 0
+dbTarget = 10000 # triggers extra output on this suite
+
+# A collection of variables used for output
+class OutNote:
+ pass
+
+outNote = OutNote()
+outNote.version = version
+outNote.comment = ""
+outNote.wannabes = 0
+outNote.outliers = 0
+
+
+# ***main()******************************************************************
+def main():
+ global dbCounter # for debugging KPB 210222
+
+ # 1. read the input
+ if args.infile != "":
+ inFile = open(args.infile)
+# elif sys.gettrace() is not None: # how to detect debugger present
+ else:
+ inFile = sys.stdin
+
+ if args.suitein or args.suitesin:
+ suites = readKinemageFile(inFile)
+ if len(suites) == 0:
+ sys.stderr.write("read no suites: perhaps wrong type of kinemage file\n")
+ sys.exit(1)
+ else:
+ residues = readResidues(inFile)
+ if len(residues) == 0:
+ sys.stderr.write("read no residues: perhaps wrong alternate code\n")
+ sys.exit(1)
+ suites = buildSuites(residues)
+ suites = suites[:-1]
+
+ # 2. process the suites
+ for s in suites:
+ if not s.validate():
+ if args.test:
+ sys.stderr.write(f"! failed validation: {s.pointID}\n")
+ suitenout.write1Suite(
+ s, bins[13], bins[13].cluster[0], 0, 0, " tangled ", "", "", "", ""
+ )
+ continue
+
+ # At this point we have a complete suite
+ bin, issue, text, pointMaster = evaluateSuite(s)
+ pointColor = "white"
+ if bin is None:
+ s.cluster = bins[0].cluster[0]
+ bins[0].cluster[0].count += 1
+ suitenout.write1Suite(
+ s, bins[0], bins[0].cluster[0], 0, 0, text, issue, "", "", ""
+ )
+ else:
+ memberPack = membership(bin, s)
+ (cluster, distance, suiteness, notes, comment,
+ pointMaster, pointColor) = memberPack
+ suitenout.write1Suite(
+ s, bin, cluster, distance, suiteness, notes, issue, comment,
+ pointMaster, pointColor)
+ s.suiteness = suiteness
+ s.distance = distance
+ s.notes = notes
+ s.pointMaster = pointMaster
+ s.pointColor = pointColor
+ dbCounter += 1
+
+ finalStats()
+ suitenout.writeFinalOutput(suites, outNote)
+
+
+# *** evaluateSuite and its tools ***************************************
+
+# Boundaries of various angle ranges
+epsilonmin = 155
+epsilonmax = 310 # 070130
+delta3min = 60
+delta3max = 105 # changed by S.J. on 06/06/2011
+delta2min = 125
+delta2max = 165
+gammapmin = 20
+gammapmax = 95 # max 070326
+gammatmin = 140
+gammatmax = 215 # max 070326
+gammammin = 260
+gammammax = 335 # max 070326
+alphamin = 25
+alphamax = 335
+betamin = 50
+betamax = 290
+zetamin = 25
+zetamax = 335
+
+# triage table for yes-no angles:
+# each of these filters, applied to a suite, will provide a
+# true or false answer as to whether this angle is in a reasonable range.
+# pointMaster is for grouping points in kinemage display
+
+# data per line: angle index, min, max, code, text, pointMaster
+triageFilters = {
+ "epsilon": (2, epsilonmin, epsilonmax, Issue.EPSILON_M, "e out", "E"),
+ "alpha": (4, alphamin, alphamax, Issue.ALPHA, " a out", "T"),
+ "beta": (5, betamin, betamax, Issue.BETA, " b out", "T"),
+ "zeta": (3, zetamin, zetamax, Issue.ZETA_M, " z out", "T"),
+}
+
+def triage(selector, suite):
+ # if angle lies outside the acceptable range, triage immediately
+ filter = triageFilters[selector]
+ index, min, max, failCode, failText, pointMaster = filter
+ if suite.angle[index] < min or suite.angle[index] > max:
+ return False, failCode, failText, pointMaster
+ else:
+ return True, None, None, ""
+
+
+# The more complex angles are handled by a "sieve".
+# A sieve will determine whether an angle is within one of several ranges
+# and provide an appropriate code indicating the range.
+# This is handled by the sift() function.
+sieveDelta = (
+ (delta3min, delta3max, 3),
+ (delta2min, delta2max, 2),
+)
+
+sieveGamma = (
+ (gammatmin, gammatmax, "t"),
+ (gammapmin, gammapmax, "p"),
+ (gammammin, gammammax, "m"),
+)
+
+def sift(sieve, angle, failCode):
+ for filter in sieve:
+ min, max, code = filter
+ if min <= angle <= max:
+ return code, "", ""
+ failMessage = failMessages[failCode]
+ return None, failCode, failMessage
+
+
+def evaluateSuite(suite):
+ global bins
+
+ # The order of triage operations, though it may seem arbitrary,
+ # was carefully chosen by the scientists.
+ ok, failCode, notes, pointMaster = triage("epsilon", suite)
+ if not ok:
+ return None, failCode, notes, pointMaster
+
+ # Angles with several meaningful ranges:
+ # for each angle, find out which range it lies in, or none
+ # this becomes a selector to help choose a bin
+ puckerdm, failCode, notes = sift(sieveDelta, suite.deltaMinus, Issue.DELTA_M)
+ if not puckerdm:
+ return None, failCode, notes, "D"
+
+ puckerd, failCode, notes = sift(sieveDelta, suite.delta, Issue.DELTA)
+ if not puckerd:
+ return None, failCode, notes, "D"
+
+ gammaname, failCode, notes = sift(sieveGamma, suite.gamma, Issue.GAMMA)
+ if not gammaname:
+ return None, failCode, notes, "T"
+
+ ok, failCode, notes, pointMaster = triage("alpha", suite)
+ if not ok:
+ return None, failCode, notes, pointMaster
+
+ ok, failCode, notes, pointMaster = triage("beta", suite)
+ if not ok:
+ return None, failCode, notes, pointMaster
+
+ ok, failCode, notes, pointMaster = triage("zeta", suite)
+ if not ok:
+ return None, failCode, notes, pointMaster
+
+ # We have pass the test: now use this information to select a bin
+ bin = bins[(puckerdm, puckerd, gammaname)]
+ # bins is an associated dictionary indexed by the triplet of three angle classifiers
+ # each unique triplet of classifiers selects one unique bin, for a total of 12 bins.
+ return bin, None, None, ""
+
+
+# ***membership()***************************************************************
+
+# cluster membership:
+# given the bin, we are looking for the correct cluster
+def membership(bin, suite):
+ matches = np.full(MAX_CLUSTERS, 999.9)
+ matchCount = 0
+ comment = ""
+ pointMaster = ""
+ pointColor = "white"
+
+ lDominant = bin.dominant > 0
+ if lDominant: # this bin has a dominant cluster, note it
+ dominantJ = bin.dominant
+ domCluster = bin.cluster[bin.dominant]
+
+ # find the closest cluster
+ # search every cluster in the bin except cluster 0, which is for outliers
+ closestD = 999
+ for j, c in enumerate(bin.cluster[1:], 1):
+ if c.status == "wannabe" and args.nowannabe:
+ continue
+ distance = hyperEllipsoidDistance(
+ suite.angle, bin.cluster[j].angle, 4, normalWidths
+ )
+ if distance < closestD:
+ closestD = distance
+ closestJ = j
+ closestCluster = c
+ matches[j] = distance
+ if distance < 1: # suite could be a member of this cluster
+ matchCount += 1
+
+ if matchCount == 1:
+ theCluster = closestCluster
+ situation = "1-only-one"
+
+ elif matchCount > 1 and not lDominant:
+ # dominant cluster is not a possible cluster
+ # just output than minimum distance match
+ theCluster = closestCluster
+ situation = "{matchCount}-None-dom"
+
+ elif matchCount > 1: # and lDominant
+ # find the closest cluster that is not the dominant cluster
+ closestNonD = 999
+ for j, c in enumerate(bin.cluster[1:], 1):
+ if c.status == "wannabe" and args.nowannabe:
+ continue
+ if matches[j] < closestNonD and c.dominance != "dom":
+ closestNonD = matches[j]
+ closestJ = j
+ theCluster = c
+
+ if theCluster.dominance == "sat":
+ # We need to distinguish carefully whether our suite
+ # is in the dominant or satellite cluster
+ theCluster, closestJ, situation = domSatDistinction(
+ suite, domCluster, theCluster, matches, matchCount
+ )
+ else:
+ if matches[dominantJ] < matches[closestJ]:
+ closestJ = dominantJ
+ theCluster = domCluster
+ situation = f"{matchCount}-not-sat"
+ else:
+ # no match, it's an outlier
+ closestJ = 0
+ theCluster = closestCluster
+ situation = f"outlier distance {closestD:.3}"
+ outNote.outliers += 1
+ pointMaster = "O"
+ pointColor = "white"
+
+ # final computation of suiteness
+ # this time we use all 7 dimensions
+ if dbCounter >= dbTarget and dbCounter <= dbTarget + 1: # KPB debug
+ print(suite.pointID)
+ print(suite.angle)
+ print(theCluster.name)
+
+ distance = hyperEllipsoidDistance(suite.angle, theCluster.angle, 7, normalWidths)
+ # this calculation can assign or deassign a cluster
+ if distance <= 1:
+ suiteness = (cos(pi * distance) + 1) / 2
+ if suiteness < 0.01:
+ suiteness = 0.01
+ else:
+ if closestJ != 0:
+ # 7D distance forces this suite to be an outlier
+ # so we deassign it here
+ closestJ = 0
+ comment = f"7D dist {theCluster.name}"
+ if theCluster.status == "wannabe":
+ comment += " wannabe"
+ theCluster = bin.cluster[0] # outlier
+ suiteness = 0
+
+ theCluster.count += 1
+ suite.cluster = theCluster
+ if theCluster.status == "wannabe" and not args.nowannabe:
+ outNote.wannabes = 1 # once set, stays set
+ pointColor = 0 # will be handled later!!
+ if args.test:
+ print(" [suite: %s %s 4Ddist== %f, 7Ddist== %f, suiteness==%f] \n" % \
+ (theCluster.name, suite.pointID[:11], closestD, distance, suiteness))
+ return theCluster, distance, suiteness, situation, comment, pointMaster, pointColor
+
+
+def domSatDistinction(suite, domCluster, satCluster, matches, matchCount):
+ # if dotproducts both positive, then inbetween
+ # p p
+ # dom/___sat and dom___\sat
+ closestCluster = satCluster
+ closestJ = satCluster.ordinal
+ dominantJ = domCluster.ordinal
+
+ # use vector properties of numpy.array to determine difference vectors
+ domToPoint = domCluster.angle - suite.angle
+ satToPoint = satCluster.angle - suite.angle
+ domToSat = domCluster.angle - satCluster.angle
+ satToDom = -domToSat
+
+ dps = narrowDotProduct(domToPoint, domToSat, 4)
+ spd = narrowDotProduct(satToPoint, satToDom, 4)
+
+ if dps > 0 and spd > 0:
+ # the trickiest case: point is between dom and sat
+ domWidths = normalWidths.copy()
+ if args.satellites:
+ satWidths = satelliteWidths.copy()
+ else:
+ satWidths = normalWidths.copy()
+ if satCluster.satelliteInfo is not None:
+ modifyWidths(domWidths, satWidths, satCluster.satelliteInfo)
+ disttodom = hyperEllipsoidDistance(suite.angle, domCluster.angle, 4, domWidths)
+ disttosat = hyperEllipsoidDistance(suite.angle, satCluster.angle, 4, satWidths)
+ if disttodom < disttosat:
+ closestJ = dominantJ
+ closestCluster = domCluster
+ situation = f"{matchCount}-BETWEEN-dom-sat({disttodom:7.3}|{disttosat:7.3})"
+ # else the satellite cluster remains the chosen cluster
+
+ else:
+ # the point is not in between
+ # just assign by closest standard distance evaluation
+ if matches[dominantJ] < matches[closestJ]:
+ closestJ = dominantJ
+ closestCluster = domCluster
+ situation = f"{matchCount}-OUTSIDE-dom-sat"
+
+ return closestCluster, closestJ, situation
+
+
+# *** Gathering some statistics *********************************************
+
+def finalStats():
+ for bin in bins.values():
+ for c in bin.cluster:
+ bin.count += c.count
+
+
+# *** The fancy math ********************************************************
+
+# This variable was experimental but we have settled on 3:
+power = 3
+
+def hyperEllipsoidDistance(suiteAngles, clusterAngles, nAngles, widthArray):
+ global dbCounter
+ if nAngles == 4:
+ workRange = range(2, 6)
+ else:
+ workRange = range(1, 8)
+
+ summation = 0
+ for k in workRange:
+ delta = abs(suiteAngles[k] - clusterAngles[k])
+ delta = delta / widthArray[k]
+ delToPower = pow(delta, power)
+ summation = summation + delToPower
+ if dbCounter >= dbTarget and nAngles > 4: # KPB debug 120221
+ sys.stderr.write("db=%3d, k=%d, del=%8.4f, delpower=%10.6f, dpower=%10.6f\n" %
+ (dbCounter, k, delta, delToPower, summation) )
+ result = pow(summation, 1 / power)
+ if dbCounter == dbTarget and nAngles > 4:
+ sys.stderr.write("final = %7.3f\n" % result)
+ return result
+
+
+# The narrow dot product involves only a subset of the dimensions,
+# either 4 or 6. In practice, only 4 is in use.
+def narrowDotProduct(a, b, nAngles):
+ if nAngles == 4:
+ return np.dot(a[2:6], b[2:6])
+ else:
+ return np.dot(a[1:8], b[1:8])
+
+
+def modifyWidths(dom, sat, satInfo):
+ for m in range(9):
+ if satInfo.satelliteWidths[m] > 0:
+ sat[m] = satInfo.satelliteWidths[m]
+ if satInfo.dominantWidths[m] > 0:
+ dom[m] = satInfo.dominantWidths[m]
+
+
+def showHelpText():
+ sys.stderr.write(
+ f"""
+Version {version}
+suitename -flags <stdin >stdout
+ or
+suitename inputfile -flags >stdout
+output flags: [ -report || -string || -kinemage ]
+default: -report
+
+input flags: [ -residuein || -suitein ]
+flags: [ -residuein [ -pointIDfields # ] ] default#=={args.pointidfields}
+ OR
+flags: [ -suitein [ -anglefields # ] ] default#=={args.anglefields}
+defaults: -residuein -pointIDfields {args.pointidfields}
+
+The -residuein format:
+label:model:chain:number:ins:type:alpha:beta:gamma:delta:epsilon:zeta
+if the file has alternate conformations, then use both -pointIDfields
+ and -altIDfield # to specify the number of pointID fields and which field
+ contains the altID
+use -altIDval <altID> to specify which alternate conformation to calculate
+ suite for. By default calculated for alt A
+
+-suitein takes a kinemage format, and uses records from
+ @ballists and/or @dotlists in this format:
+{{ptID}} [chi] deltam epsilon zeta alpha beta gamma delta [chi]
+ @dimension in the file, if present, overrides -anglefields
+
+Note dangle trick to make theta,...,eta suites directly
+
+flag: -report [ -chart ]
+ suites in order of input, suiteness summary at end
+( -chart : NO summary at end, for MolProbity multichart)
+
+flag: -string [-nosequence] [-oneline]
+ 3 character per suite string in order of input
+ 20 per line, ptID of n*20th at end of line
+ flag: -nosequence
+ only suite names, no Base sequence character
+ flag: -oneline
+ string all one line, no point IDs
+
+flag: -kinemage
+ kinemage of clusters grouped by pucker,pucker ...
+ group {{delta,delta}},subgroup {{gamma}},list {{cluster name}}
+ flag: -etatheta or -thetaeta
+ kinemage labels theta,eta instead of chi-1,chi
+
+flag: -satellite
+ use special general case satellite widths
+flag: -nowannabe
+ never assign suites to wannabe clusters
+Note: any DNA residues found in the input will be ignored.
+""")
+
+main()
=====================================
python/suitenamedefs.py
=====================================
@@ -0,0 +1,211 @@
+import suiteninit
+
+import numpy as np
+from numpy import array
+from enum import Enum
+
+
+# reasons why a suite may fail to be classified:
+Issue = Enum('Issue', 'DELTA_M EPSILON_M ZETA_M ALPHA BETA GAMMA DELTA')
+reasons = {
+ Issue.DELTA_M: "delta-1",
+ Issue.EPSILON_M: "epsilon-1",
+ Issue.ZETA_M: "zeta-1",
+ Issue.ALPHA: "alpha",
+ Issue.BETA: "beta",
+ Issue.GAMMA: "gamma",
+ Issue.DELTA: "delta"
+}
+
+failMessages = {
+ Issue.DELTA_M: "bad deltam",
+ Issue.GAMMA: "g out",
+ Issue.DELTA: "bad delta"
+}
+
+
+# primary (coarse grained) classification of suites
+class Bin:
+ # permanence properties
+ name = ""
+ ordinal = 0
+ cluster = ()
+ # a tuple of cluster objects
+ dominant = -1
+ # statistics gathered during the run
+ count = 0
+ active = False
+
+ def __init__(self, ordinal, name, clusters=()):
+ self.ordinal = ordinal
+ self.name = name
+ self.cluster = clusters
+ self.dominant = -1
+ self.active = False
+
+ for i, c in enumerate(clusters):
+ if c.dominance == "dom":
+ self.dominant = i
+ break
+
+
+# secondary (fine grained) classification of suite
+class Cluster:
+ # intrinsic data:
+ ordinal = 0 # its place in the bin
+ name = "" # the name of the cluster
+ status = "" # certain, wannabe, triaged, outlier, nothing, incomplete
+ clusterColor = "" # kinemage color names
+ dominance = "" # dom, sat, ord, out, tri, inc
+ satelliteInfo = None # present only if this cluster is a satellite
+
+ # tuple of 9 angles: chi-1 as 0 and chi as 8
+ # the standard 7 angles are indices 1-7:
+ angle = ()
+
+ # gathered statistics:
+ count = 0 # number of data points found in this cluster
+ suitenessSum = 0
+ suitenessCounts = None
+
+ def __init__(self, ordinal, name, status, color, dominance, angles):
+ self.ordinal = ordinal
+ self.name = name
+ self.LOK = (name != "!!")
+ self.status = status
+ self.clusterColor = color
+ self.dominance = dominance
+ self.angle = array(angles)
+ if self.dominance == "sat":
+ self.satelliteInfo = suiteninit.getSatelliteInfo(name)
+ else:
+ self.satelliteInfo = None
+ self.suitenessCounts = np.zeros(12)
+ self.suitenessSum = 0
+
+
+class SatelliteInfo:
+ # numbers used when suite is between satellite and dominant centers
+ name = ""
+ satelliteWidths = () # vector of 9 angles
+ dominantWidths = () # vector of 9 angles
+
+ def __init__(self, name, satelliteWidths, dominantWidths):
+ self.name = name
+ self.satelliteWidths = satelliteWidths
+ self.dominantWidths = dominantWidths
+
+
+# Residue: a raw residue as normally represented
+class Residue:
+ '''
+ A residue as normally read in.
+ Used only briefly as input.
+ '''
+ pointIDs = []
+ base = " " # A, C, G, U, ...
+ # The 6 angles:
+ alpha = 0
+ beta = 0
+ gamma = 0
+ delta = 0
+ epsilon = 0
+ zeta = 0
+ chi = 0
+ angle = np.full(7, 0.0)
+
+ def __init__(self, ID, base, angles):
+ self.pointIDs = ID
+ self.base = base
+ self.angle = angles
+ self.unpackAngles()
+
+ def unpackAngles(self):
+ self.alpha = self.angle[0]
+ self.beta = self.angle[1]
+ self.gamma = self.angle[2]
+ self.delta = self.angle[3]
+ self.epsilon = self.angle[4]
+ self.zeta = self.angle[5]
+ # for the future:
+ # we can accept chi as a seventh angle if provided
+ if len(self.angle) > 6:
+ self.chi = self.angle[6]
+ else:
+ self.chi = -431602080.00 #180
+ # A preposterous compromise with the past for now
+ #self.chi = -180
+
+
+# Suite: the
+class Suite:
+ '''
+ The set of angles forming the linkage BETWEEN residues.
+ This is the core data structure used in most operations of the program.
+ '''
+ pointID = ()
+ base = " " # A, C, G, U, ...
+ #chiMinus = 0
+ deltaMinus = 0
+ epsilon = 0
+ zeta = 0
+ alpha = 0
+ beta = 0
+ gamma = 0
+ delta = 0
+ chi = 0
+ # dual representation: individual angles are named for clarity
+ # array is for convenience of computation
+ angle = np.full(9, 0.0)
+
+ @property
+ def chiMinus(self):
+ return self.angle[0]
+
+ @chiMinus.setter
+ def chiMinus(self, value):
+ self.angle[0] = value
+
+ # fields computed during analysis:
+ cluster = None # The cluster to which it is assigned
+ suiteness = 0.0
+ distance = 0.0
+ notes = ""
+ pointMaster = ""
+ pointColor = ""
+
+ def __init__(self, ID, base, angles=None):
+ self.pointID = ID
+ self.base = base
+ if angles is not None:
+ self.angle = angles
+ self.unpackAngles()
+ self.cluster = None
+ self.suiteness = 0.0
+ self.distance = 0.0
+ self.notes = ""
+
+ def validate(self):
+ # make sure that angles deltaMinus through delta are reasonable
+ for i in range(1, 8):
+ if self.angle[i] < 0 or self.angle[i] > 360:
+ return False
+ return True
+
+ def gatherAngles(self):
+ self.angle = array((self.chiMinus, self.deltaMinus, self.epsilon, self.zeta,
+ self.alpha, self.beta, self.gamma, self.delta, self.chi))
+
+ def unpackAngles(self):
+ self.chiMinus = self.angle[0]
+ self.deltaMinus = self.angle[1]
+ self.epsilon = self.angle[2]
+ self.zeta = self.angle[3]
+ self.alpha = self.angle[4]
+ self.beta = self.angle[5]
+ self.gamma = self.angle[6]
+ self.delta = self.angle[7]
+ self.chi = self.angle[8]
+
+
+
=====================================
python/suiteninit.py
=====================================
@@ -0,0 +1,391 @@
+import suitenamedefs
+from suitenamedefs import Bin, Cluster, SatelliteInfo
+
+from numpy import array
+import argparse, sys
+
+# suiteninit.py:
+# This is a self initializing module that embodies the data around which
+# this program is built. It exports to primary data structures:
+# args the command line arguments, stored in a Namespace
+# bins the bin and cluster definitions
+
+
+MAX_CLUSTERS = 16 # practical, observed limit of clusters in a bin
+args = {}
+
+
+# ***parseCommandLine()*******************************************************
+def parseCommandLine():
+ global args
+ for i, arg in enumerate(sys.argv):
+ sys.argv[i] = arg.lower()
+
+ parser = argparse.ArgumentParser()
+ # the input file may be given as an argument or as a redirect
+ parser.add_argument("infile", nargs="?", default="")
+
+ # input styles
+ parser.add_argument("--residuein", "-residuein", action="store_true")
+ parser.add_argument("--residuesin", "-residuesin", action="store_true")
+ parser.add_argument("--suitein", "-suitein", action="store_true")
+ parser.add_argument("--suitesin", "-suitesin", action="store_true")
+
+ # output styles (default is --report)
+ outputStyle = parser.add_mutually_exclusive_group()
+ outputStyle.add_argument("--report", "-report", action="store_true")
+ outputStyle.add_argument("--string", "-string", action="store_true")
+ outputStyle.add_argument("--kinemage", "-kinemage", action="store_true")
+ parser.add_argument(
+ "--chart", "-chart", action="store_true"
+ ) # a modifier to report
+
+ # additional options
+ parser.add_argument("--satellites", "-satellites", action="store_true")
+ parser.add_argument("--nowannabe", "-nowannabe", action="store_true")
+ parser.add_argument("--nosequence", "-nosequence", action="store_true")
+ parser.add_argument("--thetaeta", "-thetaeta", action="store_true")
+ parser.add_argument("--etatheta", "-etatheta", action="store_true")
+ parser.add_argument("--test", "-test", action="store_true")
+
+ # numerical options
+ parser.add_argument("--anglefields", "-anglefields", type=int, default=9)
+ parser.add_argument("--pointidfields", "-pointidfields", type=int, default=6)
+ parser.add_argument("--ptid", "-ptid", type=int, default=0)
+ parser.add_argument("--altid", "-altid", type=str, default="A")
+ parser.add_argument("--altidval", "-altidval", type=str, default="A")
+ parser.add_argument("--altidfield", "-altidfield", type=int, default=4)
+
+ # the following are deprecated:
+ parser.add_argument("--angles", type=int, default=9)
+ parser.add_argument("--resAngles", type=int, default=6)
+
+ # now actually parse them
+ args = parser.parse_args()
+ if args.ptid:
+ args.pointidfields = args.ptid
+ pass
+
+
+# *** codes to match various residues ****************************
+
+idFields = 5
+match_list = (
+ (":ADE: A:A : Ar:ATP:ADP:AMP:T6A:1MA:RIA: I:I :", "A"),
+ (":GUA: G:G : Gr:GTP:GDP:GMP:GSP:1MG:2MG:M2G:OMG: YG: 7MG:YG :", "G"),
+ (":CYT: C:C : Cr:CTP:CDP:CMP:5MC:OMC:", "C"),
+ (":URA:URI: U: Ur:U :UTP:UDP:UMP:5MU:H2U:PSU:4SU:", "U"),
+ (":THY: T:T : Tr:TTP:TDP:TMP:", "T"),
+)
+
+
+# *** Special data relating to satellite clusters ****************
+
+# This function operates on the satelliteData list below
+# it creates an associated dictionary based on the name
+def buildSatelliteTable():
+ global satelliteTable
+ satelliteTable = {}
+ for item in satelliteData:
+ name = item[0]
+ satWidths = item[1]
+ domWidths = item[2]
+ satelliteTable[name] = SatelliteInfo(name, satWidths, domWidths)
+
+
+# The satellite data:
+# The widths below are used to determine multidimensional hyperellipsoidal distances.
+# A distance <= 1 is considered "in" the cluster
+
+# There are three tiers of widths here:
+# 1. The normal widths, deltamw etc.
+# 2. The general satellite widths, epsilonsatw etc.
+# 3. The special satellite widths in the table farther below, satelliteData
+# The normalWidths are used for a typical cluster.
+# The satelliteWidths are an exception, used for satellite clusters.
+# The satelliteData are exceptions to the satelliteWidths, for certain specific
+# satellite clusters.
+
+# SITUATION as of 210213:
+# The satelliteWidths are used ONLY if the --satellites arg is used
+
+clusterhalfwidthsversion = "070328"
+deltamw = 28
+epsilonw = 60
+epsilonsatw = 50 # satw 070328
+zetaw = 55
+zetasatw = 50 # satw 070328
+alphaw = 50
+alphasatw = 45 # satw 070328
+betaw = 70
+betasatw = 60 # satw 070328
+gammaw = 35
+deltaw = 28
+
+
+# width arrays set the widths of clusters in the various dimensions
+# the zeroes on either end may someday be replaced with widths for the
+# chi angles
+normalWidths = array((0, deltamw, epsilonw, zetaw, alphaw, betaw, gammaw, deltaw, 0))
+satelliteWidths = array(
+ (0, deltamw, epsilonsatw, zetasatw, alphasatw, betasatw, gammaw, deltaw, 0)
+)
+
+satelliteData = (
+ # sat 9 angles sat widths 9 angles dom width
+ ("1m", (0, 0, 0, 0, 0, 32, 0, 0, 0), (0, 0, 0, 0, 0, 64, 0, 0, 0)),
+ ("1L", (0, 0, 18, 0, 0, 18, 0, 0, 0), (0, 0, 70, 0, 0, 70, 0, 0, 0)),
+ ("&a", (0, 0, 20, 20, 0, 0, 0, 0, 0), (0, 0, 60, 60, 0, 0, 0, 0, 0)),
+ ("1f", (0, 0, 0, 0, 0, 47, 0, 0, 0), (0, 0, 0, 0, 0, 65, 0, 0, 0)),
+ ("1[", (0, 0, 0, 0, 0, 34, 0, 0, 0), (0, 0, 0, 0, 0, 56, 0, 0, 0)),
+ ("4a", (0, 0, 40, 40, 0, 0, 0, 0, 0), (0, 0, 50, 50, 0, 0, 0, 0, 0)),
+ ("#a", (0, 0, 26, 26, 0, 0, 0, 0, 0), (0, 0, 36, 36, 0, 0, 0, 0, 0)),
+ ("0i", (0, 0, 0, 0, 0, 60, 0, 0, 0), (0, 0, 0, 0, 0, 60, 0, 0, 0)),
+ ("6j", (0, 0, 0, 0, 0, 60, 0, 0, 0), (0, 0, 0, 0, 0, 60, 0, 0, 0)),
+)
+# note on the satelliteData:
+# The satellite cluster is stated. The dominant cluster can be found by looking
+# in the bin containing the satellite cluster, it will be the cluster with a
+# dominance of "dom"
+
+
+def getSatelliteInfo(name):
+ return satelliteTable[name]
+
+
+# *** Cluster data *******************************************************
+
+# The cluster data: centers of each cluster in 7 dimensions
+# (number, name, status, color, dominance ... the 7 angles)
+bin0data = (0, "trig",
+ ( 0 , "!!", "triaged", "white ", "tri",
+ (0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0)),
+)
+
+bin1data = (1, "33 p",
+ ( 0 , "!!", "outlier", "white ", "out",
+ (0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0)),
+ ( 1 , "1a", "certain", "yellowtint ", "dom",
+ (180.0, 81.495, 212.25, 288.831, 294.967, 173.99, 53.55, 81.035, 180.0)),
+ ( 2 , "1m", "certain", "blue ", "sat",
+ (180.0, 83.513, 218.12, 291.593, 292.247, 222.3, 58.067, 86.093, 180.0)),
+ ( 3 , "1L", "certain", "green ", "sat",
+ (180.0, 85.664, 245.014, 268.257, 303.879, 138.164, 61.95, 79.457, 180.0)),
+ ( 4 , "&a", "certain", "cyan ", "sat",
+ (180.0, 82.112, 190.682, 264.945, 295.967, 181.839, 51.455, 81.512, 180.0)),
+ ( 5 , "7a", "certain", "pink ", "ord",
+ (180.0, 83.414, 217.4, 222.006, 302.856, 160.719, 49.097, 82.444, 180.0)),
+ ( 6 , "3a", "certain", "magenta ", "ord",
+ (180.0, 85.072, 216.324, 173.276, 289.32, 164.132, 45.876, 84.956, 180.0)),
+ ( 7 , "9a", "certain", "hotpink ", "ord",
+ (180.0, 83.179, 210.347, 121.474, 288.568, 157.268, 49.347, 81.047, 180.0)),
+ ( 8 , "1g", "certain", "sea ", "ord",
+ (180.0, 80.888, 218.636, 290.735, 167.447, 159.565, 51.326, 85.213, 180.0)),
+ ( 9 , "7d", "certain", "purple ", "ord",
+ (180.0, 83.856, 238.75, 256.875, 69.562, 170.2, 52.8, 85.287, 180.0)),
+ ( 10 , "3d", "certain", "peach ", "ord",
+ (180.0, 85.295, 244.085, 203.815, 65.88, 181.13, 54.68, 86.035, 180.0)),
+ ( 11 , "5d", "certain", "yellow ", "ord",
+ (180.0, 79.671, 202.471, 63.064, 68.164, 143.45, 49.664, 82.757, 180.0)),
+ ( 12 , "3g", "wannabe", "gray ", "ord",
+ (180.0, 84.0, 195.0, 146.0, 170.0, 170.0, 52.0, 84.0, 180.0)),
+)
+
+bin2data = (2, "33 t",
+ ( 0 , "!!", "outlier", "white ", "out",
+ (0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0)),
+ ( 1 , "1e", "certain", "red ", "ord",
+ (180.0, 80.514, 200.545, 280.51, 249.314, 82.662, 167.89, 85.507, 180.0)),
+ ( 2 , "1c", "certain", "gold ", "dom",
+ (180.0, 80.223, 196.591, 291.299, 153.06, 194.379, 179.061, 83.648, 180.0)),
+ ( 3 , "1f", "certain", "lime ", "sat",
+ (180.0, 81.395, 203.03, 294.445, 172.195, 138.54, 175.565, 84.47, 180.0)),
+ ( 4 , "5j", "certain", "sky ", "ord",
+ (180.0, 87.417, 223.558, 80.175, 66.667, 109.15, 176.475, 83.833, 180.0)),
+ ( 5 , "5n", "wannabe", "gray ", "ord",
+ (180.0, 86.055, 246.502, 100.392, 73.595, 213.752, 183.395, 85.483, 180.0)),
+)
+
+bin3data = (3, "33 m",
+ ( 0 , "!!", "outlier", "white ", "out",
+ (0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0)),
+ ( 1 , "!!", "nothing", "white ", "out",
+ (0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0)),
+)
+
+bin4data = (4, "32 p",
+ ( 0 , "!!", "outlier", "white ", "out",
+ (000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000)),
+ ( 1 , "1b", "certain", "cyan ", "dom",
+ (180.000, 084.215, 215.014, 288.672, 300.420, 177.476, 058.307, 144.841, 180.000)),
+ ( 2 , "1[", "certain", "pink ", "sat",
+ (180.000, 082.731, 220.463, 288.665, 296.983, 221.654, 054.213, 143.771, 180.000)),
+ ( 3 , "3b", "certain", "lilac ", "ord",
+ (180.000, 084.700, 226.400, 168.336, 292.771, 177.629, 048.629, 147.950, 180.000)),
+ ( 4 , "1z", "certain", "peach ", "ord",
+ (180.000, 083.358, 206.042, 277.567, 195.700, 161.600, 050.750, 145.258, 180.000)),
+ ( 5 , "5z", "certain", "purple ", "ord",
+ (180.000, 082.614, 206.440, 052.524, 163.669, 148.421, 050.176, 147.590, 180.000)),
+ ( 6 , "7p", "certain", "sea ", "ord",
+ (180.000, 084.285, 236.600, 220.400, 068.300, 200.122, 053.693, 145.730, 180.000)),
+ ( 7 , "5p", "wannabe", "gray ", "ord",
+ (180.000, 084.457, 213.286, 069.086, 075.500, 156.671, 057.486, 147.686, 180.000)),
+)
+
+bin5data = (5, "32 t",
+ ( 0 , "!!", "outlier", "white ", "out",
+ (000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000)),
+ ( 1 , "1t", "certain", "red ", "ord",
+ (180.000, 081.200, 199.243, 288.986, 180.286, 194.743, 178.200, 147.386, 180.000)),
+ ( 2 , "5q", "certain", "yellow ", "ord",
+ (180.000, 082.133, 204.933, 069.483, 063.417, 115.233, 176.283, 145.733, 180.000)),
+)
+
+bin6data = (6, "32 m",
+ ( 0 , "!!", "outlier", "white ", "out",
+ (000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000)),
+ ( 1 , "1o", "certain", "sky ", "ord",
+ (180.000, 083.977, 216.508, 287.192, 297.254, 225.154, 293.738, 150.677, 180.000)),
+ ( 2 , "7r", "certain", "lilactint ", "ord",
+ (180.000, 084.606, 232.856, 248.125, 063.269, 181.975, 295.744, 149.744, 180.000)),
+ ( 3 , "5r", "wannabe", "gray ", "ord",
+ (180.000, 083.000, 196.900, 065.350, 060.150, 138.425, 292.550, 154.275, 180.000)),
+)
+
+bin7data = (7, "23 p",
+ ( 0 , "!!", "outlier", "white ", "out",
+ (000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000)),
+ ( 1 , "2a", "certain", "cyan ", "ord",
+ (180.000, 145.399, 260.339, 288.756, 288.444, 192.733, 053.097, 084.067, 180.000)),
+ ( 2 , "4a", "certain", "yellow ", "sat",
+ (180.000, 146.275, 259.783, 169.958, 298.450, 169.583, 050.908, 083.967, 180.000)),
+ ( 3 , "0a", "certain", "green ", "dom",
+ (180.000, 149.286, 223.159, 139.421, 284.559, 158.107, 047.900, 084.424, 180.000)),
+ ( 4 , "#a", "certain", "hotpink ", "sat",
+ (180.000, 148.006, 191.944, 146.231, 289.288, 150.781, 042.419, 084.956, 180.000)),
+ ( 5 , "4g", "certain", "greentint ", "ord",
+ (180.000, 148.028, 256.922, 165.194, 204.961, 165.194, 049.383, 082.983, 180.000)),
+ ( 6 , "6g", "certain", "gold ", "ord",
+ (180.000, 145.337, 262.869, 079.588, 203.863, 189.688, 058.000, 084.900, 180.000)),
+ ( 7 , "8d", "certain", "red ", "ord",
+ (180.000, 148.992, 270.596, 240.892, 062.225, 176.271, 053.600, 087.262, 180.000)),
+ ( 8 , "4d", "certain", "sky ", "ord",
+ (180.000, 149.822, 249.956, 187.678, 080.433, 198.133, 061.000, 089.378, 180.000)),
+ ( 9 , "6d", "certain", "orange ", "ord",
+ (180.000, 146.922, 241.222, 088.894, 059.344, 160.683, 052.333, 083.417, 180.000)),
+ ( 10 , "2g", "wannabe", "gray ", "ord",
+ (180.000, 141.900, 258.383, 286.517, 178.267, 165.217, 048.350, 084.783, 180.000)),
+)
+
+bin8data = (8, "23 t",
+ ( 0 , "!!", "outlier", "white ", "out",
+ (000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000)),
+ ( 1 , "2h", "certain", "sea ", "ord",
+ (180.000, 147.782, 260.712, 290.424, 296.200, 177.282, 175.594, 086.565, 180.000)),
+ ( 2 , "4n", "certain", "peach ", "ord",
+ (180.000, 143.722, 227.256, 203.789, 073.856, 216.733, 194.444, 080.911, 180.000)),
+ ( 3 , "0i", "certain", "lilactint ", "sat",
+ (180.000, 148.717, 274.683, 100.283, 080.600, 248.133, 181.817, 082.600, 180.000)),
+ ( 4 , "6n", "certain", "lilac ", "dom",
+ (180.000, 150.311, 268.383, 084.972, 063.811, 191.483, 176.644, 085.600, 180.000)),
+ ( 5 , "6j", "certain", "purple ", "sat",
+ (180.000, 141.633, 244.100, 066.056, 071.667, 122.167, 182.200, 083.622, 180.000)),
+)
+
+bin9data = (9, "23 m",
+ ( 0 , "!!", "outlier", "white ", "out",
+ (000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000)),
+ ( 1 , "0k", "wannabe", "gray ", "ord",
+ (180.000, 149.070, 249.780, 111.520, 278.370, 207.780, 287.820, 086.650, 180.000)),
+)
+
+bin10data = (10, "22 p",
+ ( 0 , "!!", "outlier", "white ", "out",
+ (000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000)),
+ ( 1 , "2[", "certain", "sea ", "ord",
+ (180.000, 146.383, 259.402, 291.275, 291.982, 210.048, 054.412, 147.760, 180.000)),
+ ( 2 , "4b", "certain", "gold ", "ord",
+ (180.000, 145.256, 244.622, 162.822, 294.159, 171.630, 045.900, 145.804, 180.000)),
+ ( 3 , "0b", "certain", "red ", "ord",
+ (180.000, 147.593, 248.421, 112.086, 274.943, 164.764, 056.843, 146.264, 180.000)),
+ ( 4 , "4p", "certain", "purple ", "ord",
+ (180.000, 150.077, 260.246, 213.785, 071.900, 207.638, 056.715, 148.131, 180.000)),
+ ( 5 , "6p", "certain", "sky ", "ord",
+ (180.000, 146.415, 257.831, 089.597, 067.923, 173.051, 055.513, 147.623, 180.000)),
+ ( 6 , "2z", "wannabe", "gray ", "ord",
+ (180.000, 142.900, 236.550, 268.800, 180.783, 185.133, 054.467, 143.350, 180.000)),
+)
+
+bin11data = (11, "22 t",
+ ( 0 , "!!", "outlier", "white ", "out",
+ (000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000)),
+ ( 1 , "4s", "certain", "lime ", "ord",
+ (180.000, 149.863, 247.562, 170.488, 277.938, 084.425, 176.413, 148.087, 180.000)),
+ ( 2 , "2u", "wannabe", "gray ", "ord",
+ (180.000, 143.940, 258.200, 298.240, 279.640, 183.680, 183.080, 145.120, 180.000)),
+)
+
+bin12data = (12, "22 m",
+ ( 0 , "!!", "outlier", "white ", "out",
+ (000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000)),
+ ( 1 , "2o", "certain", "hotpink ", "ord",
+ (180.000, 147.342, 256.475, 295.508, 287.408, 194.525, 293.725, 150.458, 180.000)),
+)
+
+bin13data = (13, "inc ",
+ ( 0 , "__", "incompl", "white ", "inc",
+ (000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000, 000.000)),
+)
+
+
+def buildBin(data):
+ ordinal = data[0]
+ name = data[1]
+ clusters = []
+ for item in data[2:]:
+ c = suitenamedefs.Cluster(*item)
+ clusters.append(c)
+ bin = Bin(ordinal, name, clusters)
+ return bin
+
+
+# The bins become an associative table that maps the bin selector information
+# directly to the bin
+# selector = (puckerdm, puckerd, gammaname)
+# bins 0 and 13 are catchbasins for outliers, they are not indexed by selectors
+def buildTheBins():
+ bins = {}
+ bins[0] = buildBin(bin0data)
+ bins[(3, 3, "p")] = buildBin(bin1data)
+ bins[(3, 3, "t")] = buildBin(bin2data)
+ bins[(3, 3, "m")] = buildBin(bin3data)
+ bins[(3, 2, "p")] = buildBin(bin4data)
+ bins[(3, 2, "t")] = buildBin(bin5data)
+ bins[(3, 2, "m")] = buildBin(bin6data)
+ bins[(2, 3, "p")] = buildBin(bin7data)
+ bins[(2, 3, "t")] = buildBin(bin8data)
+ bins[(2, 3, "m")] = buildBin(bin9data)
+ bins[(2, 2, "p")] = buildBin(bin10data)
+ bins[(2, 2, "t")] = buildBin(bin11data)
+ bins[(2, 2, "m")] = buildBin(bin12data)
+ bins[13] = buildBin(bin13data)
+
+ # build aliases so that bins can be indexed by number during output
+ bins[1] = bins[(3, 3, "p")]
+ bins[2] = bins[(3, 3, "t")]
+ bins[3] = bins[(3, 3, "m")]
+ bins[4] = bins[(3, 2, "p")]
+ bins[5] = bins[(3, 2, "t")]
+ bins[6] = bins[(3, 2, "m")]
+ bins[7] = bins[(2, 3, "p")]
+ bins[8] = bins[(2, 3, "t")]
+ bins[9] = bins[(2, 3, "m")]
+ bins[10] = bins[(2, 2, "p")]
+ bins[11] = bins[(2, 2, "t")]
+ bins[12] = bins[(2, 2, "m")]
+ return bins
+
+
+parseCommandLine()
+buildSatelliteTable()
+bins = buildTheBins()
=====================================
python/suiteninput.py
=====================================
@@ -0,0 +1,207 @@
+from suitenamedefs import Suite, Residue
+from suiteninit import args
+
+import numpy as np
+import math, sys
+
+altidfield = args.altidfield # where to find codes for alternatives
+
+# The great variety of codes that may represent each base in the input file
+NAListA = ":ADE: A:A : Ar:ATP:ADP:AMP:T6A:1MA:RIA: I:I :"
+NAListG = ":GUA: G:G : Gr:GTP:GDP:GMP:GSP:1MG:2MG:M2G:OMG: YG: 7MG:YG :"
+NAListC = ":CYT: C:C : Cr:CTP:CDP:CMP:5MC:OMC:"
+NAListU = ":URA:URI: U: Ur:U :UTP:UDP:UMP:5MU:H2U:PSU:4SU:"
+NAListT = ":THY: T:T : Tr:TTP:TDP:TMP:"
+IgnoreDNAList = ": DA: DG: DC: DT:"
+
+
+# out of the noise, determine the base
+def findBase(baseCode):
+ if len(baseCode) != 3:
+ return 'Z'
+
+ if NAListA.find(baseCode) >= 0: base='A'
+ elif NAListG.find(baseCode) >= 0: base='G'
+ elif NAListC.find(baseCode) >= 0: base='C'
+ elif NAListU.find(baseCode) >= 0: base='U'
+ elif NAListT.find(baseCode) >= 0: base='T'
+ elif IgnoreDNAList.find(baseCode) >= 0:
+ return None # we ignore DNA residues
+ else:
+ base='Y'
+ return base
+
+
+def stringToFloat(string):
+ try:
+ n = float(string)
+ except ValueError:
+ n = 9999.0 # or maybe math.nan?
+ return n
+
+
+def readResidues(inFile):
+ lines = inFile.readlines()
+ residues = []
+ for line in lines:
+ if len(line.strip()) == 0 or line[0] == '#': # blank or comment line
+ continue
+ fields = line.split(':')
+ ids = fields[:args.pointidfields]
+ print(args.pointidfields, ids, line)
+
+ baseCode = fields[args.pointidfields-1]
+ angleStrings = fields[args.pointidfields:]
+ if ids[altidfield].strip() != "" and ids[altidfield] != args.altidval:
+ continue # lines for the wrong alternative conformation are ignored
+
+ base = findBase(baseCode)
+ if not base: # ignore DNA bases
+ continue
+ angles = np.array([stringToFloat(s) for s in angleStrings])
+ for i in range(len(angles)):
+ if angles[i] < 0:
+ angles[i] += 360.0
+
+ residue = Residue(ids, base, angles)
+ residues.append(residue)
+ return residues
+
+
+def readKinemageFile(inFile):
+ """
+ We glean the following information from a kinemage file:
+ The @dimension command gives us the number of dimensions in the data
+ Anything between a @balllist command and a subsequent @ command
+ is a data line.
+ """
+ lines = inFile.readlines()
+ goodLines = []
+ place, line = findPrefixInList(lines, "@dimension")
+ if place > 0:
+ items = line.split()
+ dimension = len(items) - 1
+ else:
+ dimension = args.anglefields
+ place = 0
+ while place >= 0:
+ begin, line = findPrefixesInList(lines, "@balllist", "@dotlist", place)
+ if begin > 0:
+ end, line = findPrefixInList(lines, "@", begin+1)
+ place = end
+ if end < 0: end = len(lines)
+ goodLines += lines[begin + 1 : end]
+ else:
+ break
+ if len(goodLines) == 0:
+ goodLines = lines # assume a pure data file
+ return readKinemageSuites(goodLines, dimension)
+
+
+def readKinemageSuites(lines, dimension):
+ """Read a list of kinemage data lines to yield a suite."""
+ suites = []
+ for line in lines:
+ if len(line.strip()) == 0 or line[0] == '#': # blank or comment line
+ continue
+ # A meaningful line begins with an id string enclosed in braces
+ if line[0] == '{':
+ mark = line.find('}')
+ if mark > 0:
+ idString = line[1:mark]
+ ids = idString.split(':')
+
+ # there may be some miscellaneous markers after the id string
+ k = mark + 1
+ while k < len(line) and not line[k].isdigit():
+ k = k + 1
+ mark2 = k
+
+ # once we see a number, everything else is angles
+ angleText = line[mark2:]
+ angleStrings = angleText.split(' ')
+ angleStrings2 = angleText.split(',')
+ if len(angleStrings2) > len(angleStrings):
+ angleStrings = angleStrings2
+ angleList = [stringToFloat(s) for s in angleStrings]
+ if len(angleList) != dimension:
+ continue # wrong number of dimensions means probably not a data point
+ if dimension == 9:
+ angles = np.array(angleList)
+ else: # given only 7 angles,skipping the chi angles on the ends
+ angles = np.array([180.0] + angleList + [180.0])
+ for i in range(len(angles)):
+ if angles[i] < 0:
+ angles[i] += 360.0
+
+ suite = Suite(ids, 'X', angles)
+ suites.append(suite)
+ return suites
+
+
+def findPrefixInList(list, prefix, start=0):
+ for i, s in enumerate(list[start:]):
+ if s.startswith(prefix):
+ return i + start, s
+ return -1, None
+
+def findPrefixesInList(list, prefix1, prefix2, start=0):
+ for i, s in enumerate(list[start:]):
+ if s.startswith(prefix1) or s.startswith(prefix2):
+ return i + start, s
+ return -1, None
+
+
+def buildSuiteBetweenResidues(r1, r2):
+ suite = Suite(r2.pointIDs, r2.base)
+ suite.chiMinus = r1.chi
+ suite.deltaMinus = r1.delta
+ suite.epsilon = r1.epsilon
+ suite.zeta = r1.zeta
+ suite.alpha = r2.alpha
+ suite.beta = r2.beta
+ suite.gamma = r2.gamma
+ suite.delta = r2.delta
+ suite.chi = r2.chi
+
+ suite.gatherAngles()
+ return suite
+
+
+def buildSuiteFirst(r2):
+ suite = Suite(r2.pointIDs, r2.base)
+ suite.alpha = r2.alpha
+ suite.beta = r2.beta
+ suite.gamma = r2.gamma
+ suite.delta = r2.delta
+ suite.epsilon = 999
+ suite.zeta = 999
+ suite.chiMinus = 999
+ suite.deltaMinus = 999
+
+ suite.gatherAngles()
+ return suite
+
+
+def buildSuiteLast(r1):
+ suite = Suite((),"")
+ suite.chiMinus = r1.chi
+ suite.deltaMinus = r1.delta
+ suite.epsilon = r1.epsilon
+ suite.zeta = r1.zeta
+ suite.alpha = 999
+ suite.beta = 999
+ suite.gamma = 999
+ suite.delta = 999
+ suite.chi = 999
+
+ suite.gatherAngles()
+ return suite
+
+
+def buildSuites(residues):
+ suites = [buildSuiteFirst(residues[0])]
+ for i in range(len(residues) - 1):
+ suites.append(buildSuiteBetweenResidues(residues[i], residues[i+1]))
+ suites.append(buildSuiteLast(residues[-1]))
+ return suites
=====================================
python/suitenout.py
=====================================
@@ -0,0 +1,393 @@
+from suitenamedefs import Bin, Cluster, Issue, reasons
+from suiteninit import args, bins
+
+import sys
+from enum import Enum
+import numpy as np
+from numpy import array
+
+outFile = sys.stdout # future: might provide a way to change this
+
+reportCountAll = 0
+trigCountAll = 0
+suitenessSumAll = 0
+binSuiteCountAll = 0
+
+
+def write1Suite(suite, bin, cluster, distance, suiteness, notes,
+ issue, comment, pointMaster, pointColor):
+ if args.string:
+ string1Suite(suite, cluster)
+ elif args.kinemage:
+ kinemage1Suite(suite, bin, cluster, notes, distance, suiteness,
+ issue, comment, pointMaster, pointColor,
+ )
+ else:
+ report1Suite(suite, bin, cluster, notes, distance, suiteness, issue, comment)
+
+
+def writeFinalOutput(suites, outNote):
+ if not args.satellites:
+ outNote.comment = " all general case widths, power = 3.00"
+ else:
+ outNote.comment = " special general case satellite widths, power = 3.00"
+ if args.string:
+ pass
+ elif args.kinemage:
+ kinemageFinal(suites, outNote)
+ else:
+ reportFinal(outNote)
+
+
+def string1Suite(suite, cluster):
+ suiteCount += 1
+ if args.nosequence:
+ basestring = ":"
+ else:
+ basestring = suite.base
+ outFile.write(f"{cluster.name}{basestring}")
+
+
+def report1Suite(
+ suite, bin, cluster, notes, distance, suiteness, issue, comment
+): # ? LComment, Ltriage
+ global reportCountAll, trigCountAll, suitenessSumAll, binSuiteCountAll
+
+ # 1. write one line of output for this suite
+ reason = ""
+ if issue:
+ reason = " " + reasons[issue]
+ elif comment:
+ reason = " " + comment
+ comment = ""
+ if cluster.status == "wannabe":
+ comment = " wannabe"
+ outIDs = ":".join(suite.pointID)
+ output = (
+ f"{outIDs} {bin.name} {cluster.name} {float(suiteness):5.3f}{reason}{comment}\n"
+ )
+ outFile.write(output)
+
+ # 2. gather statistics
+ reportCountAll += 1
+
+ if bin.ordinal == 0:
+ trigCountAll += 1
+ elif bin.ordinal < 13:
+ suitenessSumAll += suiteness
+ binSuiteCountAll += 1
+
+ if cluster.ordinal == 0:
+ cluster.suitenessCounts[11] += 1
+ else:
+ cluster.suitenessSum += suiteness
+ # report in statistical baskets at intervals of 0.1:
+ # everything from 0 to 0.1 goes in bucket 1
+ # ... everything from 0.9 to 1.o goes into bucket 10
+ if suiteness == 0:
+ bucket = 0
+ else:
+ bucket = 1 + int(suiteness * 10)
+ cluster.suitenessCounts[bucket] += 1
+
+
+def reportFinal(outNote):
+ if not args.chart:
+#@deleteLine:
+ outFile.write(outNote.comment + "\n")
+ suitenessAverage(0)
+ if bins[1].cluster[1].count > 0: # Aform 1a 070325
+ suitenessAverage(1)
+ suitenessAverage(2)
+
+
+def suitenessAverage(mode):
+ """
+ Gather statistics on suiteness.
+ 12 buckets:
+ one for suiteness=0
+ one each for divisions by tenths from 0 to 1
+ one for outliers
+ """
+ bucket = np.zeros(12, dtype=int)
+ sum = 0
+ average = 0
+ allCount = 0
+ excludedCluster = None
+
+ if mode == 1:
+ # cluster 1a all by itself
+ comment = " A form (1a)"
+ cluster = bins[1].cluster[1]
+ sum = cluster.suitenessSum
+ for k in range(12):
+ bucket[k] = cluster.suitenessCounts[k]
+ else:
+ if mode == 0:
+ # all complete suites
+ startWith = 0
+ comment = "For all"
+ outFile.write(
+ "Found {} complete suites derived from {} entries\n".format(
+ binSuiteCountAll + trigCountAll, reportCountAll
+ )
+ )
+ outFile.write(
+ "{} suites were triaged, leaving {} assigned to bins\n".format(
+ trigCountAll, binSuiteCountAll
+ )
+ )
+ elif mode == 2:
+ # all complete suites except cluster 1a
+ startWith = 1
+ comment = " non-1a has"
+ excludedCluster = bins[1].cluster[1]
+ for i in range(1, 13):
+ bin = bins[i]
+ # the final bin 13 is for the pseudo-suites with incomplete
+ # angles and is ignored
+ for cluster in bin.cluster[startWith:]:
+ # cluster 0 in every bin is for outliers and is ignored
+ if cluster is excludedCluster:
+ continue # mode 2: ignore cluster 1a
+ sum += cluster.suitenessSum
+ # print(
+ # f"cluster:{cluster.name} {cluster.count:2} cluster sum: {cluster.suitenessSum} sum: {sum}\n"
+ # )
+
+ for k in range(12):
+ bucket[k] += cluster.suitenessCounts[k]
+
+ allCount = np.sum(bucket)
+ if allCount > 1:
+ average = sum / allCount
+ else:
+ average = 0
+
+ outFile.write(
+#@ "{} {} suites: average suiteness == {:5.3f}\n".format(
+ "{} {} suites: average suiteness== {:5.3f} (power==3.00)\n".format(
+ comment, allCount, average
+ )
+ )
+ if mode == 0:
+ outFile.write(f"{bucket[11]:6d} suites are outliers\n")
+ outFile.write(f"{bucket[0]:6d} suites have suiteness == 0 \n")
+ outFile.write(f"{bucket[1]:6d} suites have suiteness > 0 <.1\n")
+ for k in range(2, 10):
+ outFile.write(
+ "{:6d} suites have suiteness >=.{} <.{}\n".format(bucket[k], k - 1, k)
+ )
+ outFile.write(f"{bucket[10]:6d} suites have suiteness >=.9 \n")
+
+
+
+# ***** kinemage output format *****************************************************
+
+def kinemage1Suite(suite, bin, cluster, notes, distance, suiteness, issue, comment,
+ pointMaster, pointColor):
+ suite.pointMaster = pointMaster
+ suite.pointColor = pointColor
+ suite.notes = notes
+
+
+# static text: viewing parameters
+janesviews = """
+ at viewid {d e z}
+ at zoom 1.00
+ at zslab 200
+ at ztran 0
+ at center 197.500 172.300 178.300
+ at axischoice 2 3 4
+ at matrix
+0.07196 0.11701 -0.99052 -0.00336 0.99312 0.11707 0.99740 -0.00509 0.07186
+ at 2viewid {zag front}
+ at 2zoom 1.00
+ at 2zslab 200
+ at 2ztran 0
+ at 2center 174.091 194.887 207.768
+ at 2axischoice 4 5 7
+ at 2matrix
+0.99508 -0.00018 -0.09905 -0.00135 -0.99993 -0.01172 -0.09904 0.0118 -0.99501
+ at 3viewid {a b g}
+ at 3zoom 1.00
+ at 3zslab 200
+ at 3ztran 0
+ at 3center 175.700 189.600 64.100
+ at 3axischoice 5 6 7
+ at 3matrix
+0.99955 0.000101 0.030002 0.0002 0.99995 -0.010012 -0.030001 0.010013 0.9995
+
+"""
+
+# static text: static items in the display
+kinemageFrame = """
+ at group {frame} dominant
+ at vectorlist {frame} color= white
+P 0.000 0.000 0.000 5.000 0.000 0.000
+P 35.000 0.000 0.000 40.000 0.000 0.000
+P 80.000 0.000 0.000 160.000 0.000 0.000
+P 0.000 0.000 0.000 0.000 5.000 0.000
+P 0.000 35.000 0.000 0.000 40.000 0.000
+P 0.000 80.000 0.000 0.000 160.000 0.000
+P 0.000 0.000 0.000 0.000 0.000 5.000
+P 0.000 0.000 35.000 0.000 0.000 40.000
+P 0.000 0.000 80.000 0.000 0.000 160.000
+P 200.000 0.000 0.000 280.000 0.000 0.000
+P 320.000 0.000 0.000 360.000 0.000 0.000
+P 0.000 200.000 0.000 0.000 280.000 0.000
+P 0.000 320.000 0.000 0.000 360.000 0.000
+P 0.000 0.000 200.000 0.000 0.000 280.000
+P 0.000 0.000 320.000 0.000 0.000 360.000
+ at labellist {XYZ} color= white
+{X} 20.000 -5.000 -5.000
+{X} 380.000 -5.000 -5.000
+{Y} -5.000 20.000 -5.000
+{Y} -5.000 380.000 -5.000
+{Z} -5.000 -5.000 20.000
+{Z} -5.000 -5.000 380.000
+ at labellist {mtp} color= green
+{p} 60.000 0.000 0.000
+{t} 180.000 0.000 0.000
+{m} 300.000 0.000 0.000
+{p} 0.000 60.000 0.000
+{t} 0.000 180.000 0.000
+{m} 0.000 300.000 0.000
+{p} 0.000 0.000 60.000
+{t} 0.000 0.000 180.000
+{m} 0.000 0.000 300.000
+
+"""
+
+def kinemageFinal(suites, outNote):
+ """
+ Output the content of a kinemage file
+ The 3, 2 order may seem odd, but it is a standard
+ """
+ kinemageHeader(outNote)
+ for deltaMinus in (3, 2):
+ for delta in (3, 2):
+ binGroupOut(deltaMinus, delta, suites)
+ triaged = bins[0]
+ if triaged.count > 0:
+ outFile.write("@group {triaged} dominant dimension=9 wrap=360 select off\n")
+ binOut(triaged, suites)
+
+
+def kinemageHeader(outNote):
+ """ The invariant portion of a kinemage file """
+ outFile.write(f"@text\n {outNote.version}\n {outNote.comment}\n")
+ outFile.write("@kinemage 1\n")
+ outFile.write("@onewidth\n")
+ if args.etatheta or args.thetaeta: # 070524
+ outFile.write(
+ "@dimension {theta} {delta-1} {epsilon-1} {zeta-1} {alpha} {beta} {gamma} {delta} {eta}\n"
+ )
+ else:
+ outFile.write(
+ "@dimension {chi-1} {delta-1} {epsilon-1} {zeta-1} {alpha} {beta} {gamma} {delta} {chi}\n"
+ )
+ outFile.write(
+ "@dimminmax 0.000 360.000 0.000 360.000 0.000 360.000 0.000 360.000 0.000\
+ 360.000 0.000 360.000 0.000 360.000 0.000 360.000 0.000 360.000\n"
+ )
+ if outNote.outliers:
+ outFile.write("@pointmaster 'O' {outliers}\n")
+ if outNote.wannabes:
+ outFile.write("@master {wannabees}\n")
+ outFile.write(janesviews)
+ outFile.write(kinemageFrame)
+
+ # if(LTepsilon) {outFile.write("@pointmaster 'E' {epsilon bad}\n");}
+ # else if(LTdelta || LTdeltam)
+ # {outFile.write("@pointmaster 'D' {delta bad}\n");}
+ # else if(LTzeta || LTalpha || LTbeta || LTgamma)
+ # {outFile.write("@pointmaster 'T' {various bad}\n");}
+ # if(Loutlier) {outFile.write("@pointmaster 'O' {outliers}\n");}
+
+
+def formatAngles(angles, separator):
+ strings = ["{:7.3f}".format(a) for a in angles]
+ out = separator.join(strings)
+ return out
+
+
+def binGroupOut(deltaMinus, delta, suites):
+ # If any bin in the group has data, generate a group header
+ groupCount = 0
+ for gamma in ("p", "t", "m"):
+ bin = bins[(deltaMinus, delta, gamma)]
+ groupCount += bin.count
+ if groupCount > 0:
+ outFile.write(
+ f"@group {{{deltaMinus}{delta}}} recessiveon dimension=9"
+ " wrap=360 select animate off\n")
+
+ # generate the data
+ for gamma in ("p", "t", "m"):
+ bin = bins[(deltaMinus, delta, gamma)]
+ if bin.count > 0:
+ binOut(bin, suites)
+
+
+def binOut(bin, suites):
+ if any([c.count > 0 for c in bin.cluster]):
+ outFile.write(f"@subgroup {{{bin.name}}} recessiveon \n")
+ for cluster in bin.cluster[1:]:
+ # the first cluster, for outliers, will be handled later
+ if cluster.count > 0:
+ extras = ""
+ if cluster.status == "wannabe":
+ extras = " master= {wannabees}"
+ # display a ball for each for each suite in this cluster
+ ballList = (
+ "@balllist {{{} {}}} color= {} radius= 1 "
+ "nohilite master= {{data}}{}\n"
+ ).format(bin.name, cluster.name, cluster.clusterColor, extras)
+ # display a ring surrounding the center of the cluster
+ ringList = (
+ "@ringlist {{{} {}}} color= {} radius= 10 width= 1 "
+ "nobutton master= {{avsigma}}{}\n"
+ ).format(bin.name, cluster.name, cluster.clusterColor, extras)
+ angleList = formatAngles(cluster.angle[1:-1], ' ')
+ ringList2 = "{{{} {}}} 180 {} 180\n".format(
+ bin.name, cluster.name, angleList
+ )
+ labelList = (
+ "@labellist {{{} {}}} color= {} nobutton "
+ "master= {{labels}}{}\n"
+ ).format(bin.name, cluster.name, cluster.clusterColor, extras)
+
+ outFile.write(ballList)
+ outPoints(bin, cluster, suites, "")
+ outFile.write(ringList)
+ outFile.write(ringList2)
+ outFile.write(labelList)
+ outFile.write(ringList2)
+
+ # handle outliers if there are any:
+ if bin.cluster[0].count > 0:
+ cluster = bin.cluster[0]
+ ballList = (
+ "@balllist {{{} {}}} color= {} radius= 1 "
+ "nohilite master= {{data}}\n"
+ ).format(bin.name, cluster.name, cluster.clusterColor)
+ outFile.write(ballList)
+ outPoints(bin, cluster, suites, "'O' white")
+#@ outPoints(bin, cluster, suites, "")
+
+
+def outPoints(bin, cluster, suites, extra1):
+ for s in suites:
+ if s.cluster is cluster:
+ if bin.name == "trig": # the triage bin is especially handled
+ extra = f"'{s.pointMaster}'"
+ else:
+ extra = extra1
+ ids = ':'.join(s.pointID)
+ line = \
+ (f"{{{bin.name} {cluster.name} :D=={s.distance:5.3f}"
+ f":S=={s.suiteness:5.3f}: {ids}}} {extra} ,") \
+ + formatAngles(s.angle, ",") + "\n"
+ outFile.write(line)
+
=====================================
python/suitenutil.py
=====================================
@@ -0,0 +1,33 @@
+# This file is no longer in use, until further notice
+
+
+# This variable has been experimental but we have settled on 3
+power = 3
+
+
+def hyperEllipsoidDistance(suiteAngles, clusterAngles, nAngles, widthArray):
+ if nAngles == 4: workRange = range(2,6)
+ else: workRange = range(1,8)
+
+ summation = 0
+ for k in workRange:
+ delta = abs(suiteAngles[k] - clusterAngles[k])
+ delta = delta/widthArray[k]
+ delToPower = pow(delta,power)
+ summation = summation+delToPower
+ result=pow(summation, 1/power)
+ return result
+
+
+def hyperEllipsoidDistance2(suiteAngles, clusterAngles, nAngles, widthArray):
+ if nAngles == 4: workRange = range(2,6)
+ else: workRange = range(1,8)
+
+ del2Power = [pow(abs(suiteAngles[k] - clusterAngles[k])/widthArray[k], power) for k in workRange]
+ summation = sum(del2Power)
+ result=pow(summation, 1/power)
+ return result
+
+
+
+
View it on GitLab: https://salsa.debian.org/med-team/suitename/-/compare/a090327f09e23ee61db7b18a2607d28565681bef...e3d67f0849270539abae50c46bbfe1d4eb09eee3
--
View it on GitLab: https://salsa.debian.org/med-team/suitename/-/compare/a090327f09e23ee61db7b18a2607d28565681bef...e3d67f0849270539abae50c46bbfe1d4eb09eee3
You're receiving this email because of your account on salsa.debian.org.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://alioth-lists.debian.net/pipermail/debian-med-commit/attachments/20220128/4d3841c5/attachment-0001.htm>
More information about the debian-med-commit
mailing list