diff --git a/README b/README new file mode 100644 index 0000000..31b6ce7 --- /dev/null +++ b/README @@ -0,0 +1,110 @@ += dxTuber = + +dxTuber is the first cavity detection tool that takes into account protein and cavity dynamics +using solvent and protein residence probabilities derived from molecular dynamics (MD) simulations. +dxTuber supports all MD packages that VMD can process. + +dxTuber results can be exported in PDB format where each cavity is an individually selectable object +of coherent voxels written as pseudo atoms, each holding the residence probability information, +in the form of mass-weighted solvent densities encoded as formal b-factors. + +It calculates pore and tunnel profiles along principal axis, +furthermore cavity volumes and solvent densities can be directly extracted from the PDB files. +Based on this data cavity statistics can be obtained easily. + +dxTuber is licensed under the GPLv2 (-> gpl2-2.0.txt) + + += Dependencies = + +- GTK+ / Bonobo +- Python 2.5 and 2.6 are tested (please mail me if you tested other versions successfull ;) +- Python-gtk + - Required for GUI version +- Python-argparse for python < 2.7 (The argparse lib should be available in python 2.7+) + http://code.google.com/p/argparse/ + - Required for cmd version + + += Installation = + +1. extract and copy all files in a directory you like + (e.g. /opt/dxTuber) +2. create in your home folder a .dxTuber folder (~/.dxTuber) +3. copy dxTuber.conf into ~/.dxTuber/ +4. modify dxTuber.conf (full path of [1.] e.g. /opt/dxTuber) +5. create softlinks of + dxTuber_gtk.py | gui version + dxTuber_cmd.py | command line version + analyzeCavity.py | cmd version of analyze function + dx2pdb.py | OpenDX -> PDB convertor + pdb2dx.py | PDB -> OpenDX convertor // supports only previously processed dxTuber PDB files. + postgroupcavities.py | Pythonscript that devides cavities into subcavities + in /usr/bin (optional) +6. start dxTuber_gtk.py + + += How to use dxTuber = + +The program is started from the console by typing dxTuber_gtk.py. +While the dxTuber is driven by a graphical user interface, general status messages +and information on calculation progress are printed to the console. Protein and water +OpenDX VolMap files are imported via the OpenDX panel and specified in the Files tab. +Scan initiates a cavity search using the protein density threshold and scanning method +specified in the Settings tab. + +Once a search is completed, the found set of cavity voxels appears as a new row in the +Files tab that can be selected for clustering into coherent cavities (Group DX), +filtering (Filter), analysis (Analyze) or export. The latter can be done via the Save to DX or +Save to PDB dialogue to save a selected set of cavity voxels in PDB or OpenDX format at any time. + +Coherent cavities can be analyzed in respect to their cross-sectional area along one of +the principal axes. Therefore a principal axis and the cavity of interest need to be +specified which is done based on each cavity's unique denomination that is encoded +as atom name in the exported PDB file. + + += Tested on = + +Suse Linux 10.3 Enterprise / 11 (32 Bit) +openSuse 11.3 (32 Bit) +Ubuntu 10.4 LTS (32 Bit) + + += Known Bugs / Issues = + +- Compiz can disturb the dxTuber widget so that the main window title is not shown. +-> Disable Compiz + +- dxTuber stops reading / calculating / saving / grouping ... +-> Keep the shell in focus + +- To stop or shut down dxTuber_gtk hit 'ctrl +c' in the shell where dxTuber was startet from. +(I tried several hours to catch the term signal from gtk window with no success...) + +- After reading / scanning / saving all radio buttons are resetted +-> Ensure your selections are correct until this bug is fixed. + +- pyGTK gets outdated & can not be compiled / installed +-> use the command line version of dxTuber + + += Acknowledgements = + +### Research lab: Computational Structural Biology @ University of Bonn: ### + +Christian Kandt Supervising, Bugtesting +Nadine Fischer Bugtesting +Thomas H. Schmidt Bugtesting, Typos ;) + +### Extranal: +Mattia Sturlese Bugreport + + += Developed by = + +Martin Raunest +Contact +e-mail: raunest@bit.uni-bonn.de or m.raunest@gmx.de +xmpp: raunest@jabber.org + diff --git a/analyseDialog.libglade b/analyseDialog.libglade new file mode 100644 index 0000000..259ce09 --- /dev/null +++ b/analyseDialog.libglade @@ -0,0 +1,147 @@ + + + + + + True + 5 + 2 + + + True + True + Diameter + 0 + True + True + analyseAreaRadio + + + 1 + 2 + 3 + 4 + GTK_EXPAND + + + + + True + True + Area + 0 + True + True + + + + 1 + 2 + 2 + 3 + GTK_EXPAND + + + + + + + + + + + True + 0.039999999105930328 + Axis + GTK_JUSTIFY_RIGHT + PANGO_WRAP_WORD_CHAR + PANGO_ELLIPSIZE_START + + + GTK_FILL + + + + + + True + 0.039999999105930328 + Group + + + 1 + 2 + + + + + True + True + X-Axis + 0 + True + True + + + + 1 + 2 + GTK_EXPAND + + + + + True + True + Y-Axis + 0 + True + True + xAxisRadioButton + + + + 2 + 3 + GTK_EXPAND + + + + + True + True + Z-Axis + 0 + True + True + xAxisRadioButton + + + + 3 + 4 + GTK_EXPAND + + + + + True + True + True + Analyse + 0 + + + + 1 + 2 + 4 + 5 + GTK_EXPAND + + + + + + + diff --git a/analyzeCavity.py b/analyzeCavity.py new file mode 100755 index 0000000..fed8711 --- /dev/null +++ b/analyzeCavity.py @@ -0,0 +1,131 @@ +#!/usr/bin/env python +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +### +# +# +# v0.25 +# +# -bug in cmd parameter check block (thx for the bug report: Mattia Sturlese) +# +# +# v0.12 +# +FilterByMinDens +# +# +# v0.1 +# +cmd implementation of analyze tube +# +# +# + +import sys +import os +import re +import subprocess +from lib.analysetube import * +from lib.readpdb import * +from lib.filterdx import * + + + + +dxTuberConfFile = os.environ['HOME']+"/.dxTuber/dxTuber.conf" +conf = open (dxTuberConfFile,'r') +confLines = conf.readlines() +for line in confLines: + if re.match ('^#',line): + continue + if re.match('VERSION',line): + tmpArray = re.split('\s',line,1) + tmp = re.sub('\"','',tmpArray[1]) + tmp = re.sub('\n','',tmp) + version = re.sub('^\s','',tmp) + +#### +# Deal with arguments +#### + +###################################################################### +# If argparse is not installed global, get a static lib # +# and import it manually (static) # +# # +# http://packages.ubuntu.com/lucid/python-argparse # +# http://code.google.com/p/argparse/ # +# # +# maybe this helps for implementation ... # +# # +# sys.path.append('path_to_argparse_lib') +# from argparse import * +# # +###################################################################### + + +import argparse + +parser = argparse.ArgumentParser(description='This tool is part of the dxTuber package. It is able to analyze previously grouped cavities along principle axis x, y and z. Based on the choosen -c cavity its cross sectioneal area is plotted against the -a axis. Output will be stored in -o OUT. (1st coloumn = axis, 2nd column = area)',epilog='dxTuber v' + version) +parser.add_argument('-p','--pdb', default=False, help='dxTuber grouped outfile') +parser.add_argument('-a','--axis', default=False,choices=('x','y','z'),help='Select an axis') +parser.add_argument('-c','--cavity',default=False, help='Cavity ID (Name column in PDB file)') +parser.add_argument('-o','--out',default=False, help='Output (1st column axis 2nd column area)') +#parser.add_argument('--mindens', help='Minimum density') + + +args = parser.parse_args() +''' +args.pdb grouped outfile +args.axis axis +args.cavity cavity ID +args.out outfile +#args.mindens + +''' + +if not args.pdb or not args.axis or not args.cavity or not args.out: + sysRun = subprocess.Popen([sys.argv[0], '-h']) + sysRun.communicate() + sys.exit() + + +readPDB = ReadPDB (args.pdb) +readPDB.readPDB() +tuberObj = readPDB.getTubeDxObject() + + +if not tuberObj.getGroupes(): + sysRun = subprocess.Popen([sys.argv[0], '-h']) + sysRun.communicate() + print "Please choose a previously grouped PDB file" + print "Ensure that your cavity id start with '0'" + sys.exit() + + + +#if args.mindens: +# filterThis = FilterDx(dxObjectTube = tuberObj.getDxObject()) +# filterThis.filterByDensity( tubeDxObj= tuberObj, densityThreshold=args.mindens) +# tuberObj.setDxObject(filterThis.getDxObject()) + +analyse = AnalyseTube (tuberObj.getVoxelGroupGrid(),int(args.cavity)) +analyse.analyseTubeArea(args.axis ) +analyse.writeStatistics(args.out) + + + + + + diff --git a/changelog b/changelog new file mode 100644 index 0000000..ce5d555 --- /dev/null +++ b/changelog @@ -0,0 +1,323 @@ +Changelog + +v0.27 7. Jan 2013 + +#readPDB +-Added error message for unsupported pdb files. Example header is now printed to CLI +-bugfixes ... + +------------ + +v0.26 14. Dec 2012 + +#dxTuber_cmd +-bug /typo for 1D z scanning + +------------- + + +v0.25 22. Oct. 2012 + +#analysCavity (cmd) +-bug in parameter check block + +------------- + + +v0.24 + +#analyzecavity +-bug: pore profiles start and end now with values of 0 +-multiple zeros in one cavity are now supported AND listet in output files + usefull if disconnected cavites are in forus of a pore profile analysis + + +------------- + +v0.23 20. Jun 2012 + +#postgroup ++'core_cavities_after_filter.pdb' stores the grouped core cavities ++distcutoff parameter can now be set via cmd + +------------- + +v0.22 21. Mar 2012 + +#postgroup ++user can now select either to postgroup ALL cavites or a single cavity + (If '--cid' statement is not given) +-bugfixes + +-------------- + +v0.2 21. Feb 2012 + +#General ++pdb2dx: Converts previously previously via dxTuber created pdb files back into OpenDX file format. Usefull if cavites shall be represented as 'surface' in VMD + +#dxTuber_cmd ++supports now pdb(grid only) files as input ++fill after grouping (Experimental) ++minimum distance keep (mdk) applies an minimum distance filter and stores filtered ISV's separately for grouping // dirty dirty dirty ... ++overlap option, switches if a ISV will be stored during scanning if a protein > protMinDens is located at the very same coordinates. + +#scanning / planesearch ++'--mdk' minimum distance keep, new idea to separate cavites inside water water clusters. +-storeOverLappingISV = False <--- switches if ISV will be deleted if a protein Voxels > protMinDens is located at the very same position + +#writepdb +-ungrouped cavities now have the (atom)id '0' instead of 'DEN' -> ungrouped cavities can be processed via analyzeCavity.py | dxTuber_gui +-99999'er atom counter bug m( // for some unknown reason code was commented out... "gremlins"? +-densities up to 999999 can now be stored in pdb files + +#filter ++filterByDens: filters by voxel density values ++deleteRectangleInside: deletes voxles inside the defined rectangle //Experimental not yet testet + +#readPDB +-bugfixes, lots of... do not parse pdb files by column in future *GRML* + +#postgroup ++detected cavities can now be devided into subcavities with postgroupcavities.py + + +---------------- + +v0.11 12. August 2011 + +#General ++analyzeCavity.py + -cmd version of the analyze cavity function + +---------------- + +v0.1 10. August 2011 + +#General ++dxTuber now supports solvent files with smaler dimension than protein input files +-new dependency python-argparse http://code.google.com/p/argparse/ (in standard libs since python 2.7) + +#Bugfix +-dx2pdb files can now be imported in dxTuber_gtk + +#dxTuber_cmd ++Full parameter support +-completely rewritten +-uses python-argparse (argument parser) + +---------------- + +v0.09 5. July 2011 + +#General ++dxTuber_cmd <- first command line version of dxTuber + -uses standard parameters + -0.005 prot min dens + -0.01 water min dens + -2D scanning ++version is now stored in ~/.dxTuber/dxTuber.conf + +#writepdb ++grouped yes / no entry in PDB output header + +#planesearch +-fixed bug if constructor is called using only necessary paramemters +-minDiameter = 0 by default + +--------------- + +v0.085 17. Mar 2011 + +#General +-default density water 0.01 -> 0 +-minimum diameter for cavities during scanning 1A -> 0A + +#writepdb ++'CRYST1' entry to enable PBC view in VMD and other molecular viewer + +#Readme ++knwon bug's + +-------------- + +v0.0841 14.Mar.2011 + ++'how to use dxTuber' section in README + +-------------- + +v0.084 4.Mar.2011 + +#GUI +-changed some labels (typos) + +----------------- + +v0.083 11.Jan.2011 + +#general ++Gridsizes !=1 are now supported + +#dx2pdb ++filter by residence probabilities stored in minSolventDensity ++settings are now stored in pdb files + +#pdbfiles ++reference OpenDX filename is now saved in 'OpenDXFile' (in pdb files) + +#README ++tested on Ubuntu Lucid 10.4 LTS (32Bit) + +-------------------------- + +v0.082 22.Sep.2010 + +#readpdb ++dxTuber now supports previously processed pdb files as full compatible input files + -only files created from version 0.082+ are supported + -all settings including groupes are imported and can be analyzed further + +#GUI ++new fields in file tab + +Min Diamemter + +Min Protein density + +Min Solvent density +-made the code little bit 'nicer' + +#save to pdb ++dxTuber settings are now stored in pdb files + +version + +ProteinFile (full path) + +SolventFile (full path) + +StepSize + +ScanType + +MinDiameter + +MinProtDens + +MinSolvDens + +Applied Filter (all applied; seperated by ', ') + +#filtertube +-minor print to commandline changes + +#tubeDx ++tubeDx class has moved into lib/tubedx.py ++lots of new variables + +#README ++minor changes + +----------------- + +v0.08 16.Sep.2010 + +#new features ++User can adjust solvent minimum density via settings tab ++User can define minimum cavity diameter along priciple axis (1D,2D,3D) [Angstrom] + +#GUI-Settings tab ++Solvent density ++Minimum diameter of the cavity + +#Scanning ++additional information during scanning is given + +Protein Threshold + +Solvent Threshold + +#writepdb.py +-changed the hard coded solvent density + +#README +-minor changes + +--------------- + +v0.073 25.Aug.2010 + +#general ++README file + +Dependencies + +How to install dxTuber +#analyseDx +-all cavity profiles start at (choosen axis = 1) + +-------------- + +v0.072 20.May.2010 + +#general ++New filter 'Filter by Distance' implemented +-Fixed warnings in analysePopupWindow + - removed + - removed + - removed + +#GUI-Settings tab +-Filter by plane -> Filter by Distance ++SpinBox 'minimum Distance' in [A] + +#filterDx +-Filter by Distance: 2D filtering modified to 3D filtering + +-------------- + +v0.071 17.May.2010 + +#general +-Modified script description header (hasn't been modified since dxTuber v0.001!) +-Fixed warnings in startup of dxTuber, some deprecated stuff from glade3 were modified in dxTuber.libglade + - removed + - removed + - removed + -modified neighborspinbox initialization + -10 1 25 1 10 10 -> 10 1 25 1 10 0 + +#GUI +-For 1D scanning all xyz axis are active per default now. +-Renamed two buttons + -Analyse Tube -> Analyze + -Filter Tuber -> Filter + +--------------- + +v0.07 14.Apr.2010 + +#new features ++fill tubes after grouping (2D & 3D) + +#general +-tubeVoxelGroupGrid [x][y][z] = group1 group2 group3 + -this effects analyseTube fillTube. + -More than one group is stored here seperated by blanks as string + -Overlapping cavities are possible and still selectable via name (PDB) + +#scanning +-user gets a print on console which scan method is currently running + +#GUI ++Added checkbox for fill Tubes + +#CMD ++Print scanning method while scanning + +----------- + +v0.067 3.Mar.2010 + +#GUI +-after opening/saving/scanning 'Files'-Tab will be on focus +-after analysing cavities along princeple axis 'Analyse Tube' window will close + +#Scantube +-notification when scanning is "Done" + +#grouptube +-print progress while grouping, in percentage + +#writedx +-notification when saving is "Done" + +#writepdb ++added "TER" + "END" in tail of PDB files + + diff --git a/dx2pdb.py b/dx2pdb.py new file mode 100755 index 0000000..b5181f1 --- /dev/null +++ b/dx2pdb.py @@ -0,0 +1,89 @@ +#!/usr/bin/env python +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +### +# This script converts OpenDX files into PDB files, optionally a threshold for minimim atoms per voxel can be set. +# +# v0.084 +# +settings are stored now in pdb files +# +filter by density + +import sys, os +from lib.readopendx import * +from lib.writepdb import * +from lib.tubedx import * +from lib.buildopendx import * + + + + + +dxTuberConfFile = os.environ['HOME']+"/.dxTuber/dxTuber.conf" +conf = open (dxTuberConfFile,'r') +confLines = conf.readlines() +for line in confLines: + if re.match ('^#',line): + continue + if re.match('VERSION',line): +# print line + tmpArray = re.split('\s',line,1) + tmp = re.sub('\"','',tmpArray[1]) + tmp = re.sub('\n','',tmp) + version = re.sub('^\s','',tmp) +print "dxTuber v" + version + + +try: + sys.argv[1] + sys.argv[2] +except IndexError: + print "\nThis script converts OpenDX files into PDB files, optionally a threshold for minimim atoms per voxel can be set.\n\nusage:\ndx2pdb (minimum-residence-probabilities)\n<> Mandatory parameters \n() Optional parameters\n\n" + sys.exit() + +#print "dxTuber v" + version +print "Reading..." +readDx = ReadOpenDx (sys.argv[1]) +readDx.importDx() +#tubeDx = TubeDx (openDxObject = readDx.getDxObject(), filename = sys.argv[1],version = version,scanned = 'converted') +filter = 1 +try: + sys.argv[3] +except: + filter = 0 + +if filter == 1: + tubeDx = TubeDx (openDxObject = readDx.getDxObject(), filename = sys.argv[1],version = version,scanned = 'converted', solventThreshold = sys.argv[3]) + print "Minimum residence probability: "+sys.argv[3] + print "Start filtering..." + densities = tubeDx.getDxObject().getDensity() + densities_filtered = [[[-1 for zR in range(len (densities[0][0]))] for yR in range(len (densities[0]))] for xR in range(len (densities))] + for xCount in range (len (densities)): + for yCount in range (len (densities[0])): + for zCount in range (len (densities[0][0])): + if float(densities[xCount][yCount][zCount]) > float(sys.argv[3]): + densities_filtered[xCount][yCount][zCount] = densities[xCount][yCount][zCount] + outputLine = int ((xCount*100)/len(densities)) + print str(outputLine)+" %\r", + tmp = (BuildOpenDx (density = densities_filtered, stepSize = tubeDx.getDxObject().getStepsize(), origin = tubeDx.getDxObject().getOrigin(),dimention = tubeDx.getDxObject().getDimention())) + tubeDx.setDxObject(tmp) + print "Done" +else: + tubeDx = TubeDx (openDxObject = readDx.getDxObject(), filename = sys.argv[1],version = version,scanned = 'converted') +print "Saving..." +savePDB = WritePDB (tubeDx) +savePDB.write(sys.argv[2]) +print "Done" diff --git a/dxTuber.conf b/dxTuber.conf new file mode 100644 index 0000000..1c4f378 --- /dev/null +++ b/dxTuber.conf @@ -0,0 +1,7 @@ +###################### +# dxTuber Config file# +###################### +# +# INSTALL_DIR "location" directory where dxtuber is installed +INSTALL_DIR "/home/bit/raunest/programms/python/dxtuber" +VERSION "0.27" diff --git a/dxTuber.libglade b/dxTuber.libglade new file mode 100644 index 0000000..208dc30 --- /dev/null +++ b/dxTuber.libglade @@ -0,0 +1,721 @@ + + + + + + + + + + True + popup-menu + + + True + + + False + + + True + + + _File + True + + + True + + + gtk-new + True + True + True + + + + + gtk-open + True + True + True + + + + + + gtk-save + True + True + True + + + + + + gtk-save-as + True + True + True + + + + + True + + + + + gtk-quit + True + True + True + + + + + + + + + _Edit + True + + + True + + + gtk-cut + True + True + True + + + + + gtk-copy + True + True + True + + + + + gtk-paste + True + True + True + + + + + gtk-delete + True + True + True + + + + + + + + + _View + True + + + + + _Help + True + + + True + + + gtk-about + True + True + True + + + + + + + + + + + BONOBO_DOCK_ITEM_BEH_EXCLUSIVE | BONOBO_DOCK_ITEM_BEH_NEVER_VERTICAL | BONOBO_DOCK_ITEM_BEH_LOCKED + + + + + True + + + True + + + Open DX + True + True + True + none + + + + False + 0 + + + + + Open PDB (gridonly) + True + True + True + none + + + + False + 1 + + + + + Save to PDB + True + True + True + none + + + + False + 2 + + + + + Save to DX + True + True + True + none + + + + False + 3 + + + + + Scan + True + True + True + none + + + False + 4 + + + + + Group DX + True + True + True + none + 0.49000000953674316 + + + + False + 5 + + + + + Filter + True + True + True + none + + + + False + 6 + + + + + Analyze + True + True + True + none + 0 + + + + False + 7 + + + + + False + 0 + + + + + True + True + True + 0 + 1 + + + True + 4 + + + True + File + + + GTK_EXPAND + + + + + + True + Set as protein + + + 1 + 2 + GTK_EXPAND + + + + + + True + Set as water + + + 2 + 3 + GTK_EXPAND + + + + + + Calc avarage + + + 3 + 4 + GTK_EXPAND + + + + + + False + + + + + True + Files + + + False + tab + + + + + True + + + True + 2 + 0 + + + True + 12 + + + True + 7 + 2 + + + True + True + 0.005 + + + 1 + 2 + 3 + 4 + GTK_FILL + GTK_FILL + + + + + 1D + True + True + False + True + True + + + + GTK_FILL + + + + + 2D + True + True + False + True + True + search1DRadio + + + + 1 + 2 + GTK_FILL + + + + + 3D + True + True + False + True + True + search1DRadio + + + + 2 + 3 + GTK_FILL + + + + + Scan X + True + False + True + False + True + True + + + 1 + 2 + GTK_FILL + + + + + Scan Y + True + False + True + False + True + True + + + 1 + 2 + 1 + 2 + GTK_FILL + + + + + Scan Z + True + False + True + False + True + True + + + 1 + 2 + 2 + 3 + GTK_FILL + + + + + True + 0 + Protein density threshold + True + + + 3 + 4 + GTK_FILL + + + + + + True + 0 + Solvent density threshold + True + + + 4 + 5 + GTK_FILL + + + + + + True + True + ● + 0 + + + 1 + 2 + 4 + 5 + GTK_FILL + GTK_FILL + + + + + Fill cavities after grouping (2D & 3D) + True + True + False + True + + + 1 + 2 + 6 + 7 + GTK_FILL + + + + + True + 0 + Minimum diameter in A + True + + + 5 + 6 + GTK_FILL + + + + + + True + True + ● + 0 + + + 1 + 2 + 5 + 6 + GTK_FILL + GTK_FILL + + + + + + + + + + + + True + <b>Scan</b> + True + + + label_item + + + + + 0 + + + + + True + 0 + + + True + 12 + + + True + 2 + 2 + + + True + True + 2 0 100 1 10 0 + + + 1 + 2 + 1 + 2 + + + + + True + True + 10 1 25 1 10 0 + + + 1 + 2 + GTK_FILL + + + + + Filter by distance + True + True + False + True + True + filterNeighbourRadio + + + 1 + 2 + GTK_FILL + + + + + Neighbor + True + True + False + True + True + + + GTK_FILL + + + + + + + + + True + <b>Filter</b> + True + + + label_item + + + + + 1 + + + + + 1 + + + + + True + Settings + + + 1 + False + tab + + + + + 1 + + + + + + + True + True + + + + + True + 4 + True + True + + + 1 + + + + diff --git a/dxTuber_cmd.py b/dxTuber_cmd.py new file mode 100755 index 0000000..06717da --- /dev/null +++ b/dxTuber_cmd.py @@ -0,0 +1,471 @@ +#!/usr/bin/env python +# +# This file is part of the dxTuber package +# dxTuber v0.27 7. Jan 2013 +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +# +#################################################################################### +### This tool s written to detect cavities based on protein and solvent movement ### +### ### +### ### +### 1. flood your system with water (by using GROMACS) ### +### 2. Open your molecule in VMD ### +### 3. save 2 OpenDX (only solvent and only protein) files via VMD -> VolMap ### +### 4. start this script ### +### ### +### Contact: raunest@bit.uni-bonn.de or m.raunest@gmx.de Martin Raunest ### +### ### +#################################################################################### +# +# v.26 +# -bugfix / typo fixed fox z-axis scanning +# +# v.2 +# +fill after grouping implemented +# +pdb files as input +# +if input was not given completely or incorrect "dxtuber_cmd -h" (sys.argv[0] wuha fancy fetature ) will be printed to cmd +# + '--mdk' (minimum distance keep ... very good english i know ... ) option to separate cavities in a "solvent network" experimental ... +# +overlap, switches if a ISV will be stored during scanning if a protein > protMinDens is located at the very same coordinates. +# +# v0.1 +# -completely rewritten +# +dxTuber_cmd comes with full parameter support +# +comments :] +# +# +# v0.09 +# +wuhu dxTuber is now also on command line ... beware only standard parameters and 2D scanning are used +# +# +# + +import sys +import os +import re +import subprocess +from lib.readopendx import * +from lib.planesearch import * +from lib.writepdb import * +from lib.writedx import * +from lib.grouptubes import * +from lib.filterdx import * +from lib.analysetube import * +from lib.filltube import * +from lib.tubedx import * +from lib.readpdb import * + +######### +# Get dxTuber settings +######### + + + +dxTuberConfFile = os.environ['HOME']+"/.dxTuber/dxTuber.conf" +conf = open (dxTuberConfFile,'r') +confLines = conf.readlines() +for line in confLines: + if re.match ('^#',line): + continue + if re.match('VERSION',line): + tmpArray = re.split('\s',line,1) + tmp = re.sub('\"','',tmpArray[1]) + tmp = re.sub('\n','',tmp) + version = re.sub('^\s','',tmp) + +#### +# Deal with arguments +#### + +###################################################################### +# If argparse is not installed global, get a static lib # +# and import it manually (static) # +# # +# http://packages.ubuntu.com/lucid/python-argparse # +# http://code.google.com/p/argparse/ # +# # +# maybe this helps for implementation ... # +# # +# sys.path.append('path_to_argparse_lib') +# from argparse import * +# # +###################################################################### + + +import argparse + +parser = argparse.ArgumentParser(description='Command line version of dxTuber',epilog='dxTuber v' + version) + +### +# for static libs use: +# parser = ArgumentParser(description='Command line version of dxTuber',epilog='dxTuber v' + version) +### +parser.add_argument('-p','--protein',default=False, dest='protein',help='Protein dx file') +parser.add_argument('-s','--solvent',default=False, help='Solvent dx file') +parser.add_argument('--ipdb', default=False,help='PDB inputfile (only previously via dxTuber created pdb files are supported)') + +parser.add_argument('--sm', default='2D', choices=('1D','2D','3D'), help='Scan method (default=2D)') + +parser.add_argument('--sm_x', default=False, action='store_const', const=True, help='1D scanning x (default=False)') +parser.add_argument('--sm_y', default=False, action='store_const',const=True, help='1D scanning y (default=False)') +parser.add_argument('--sm_z', default=False, action='store_const', const=True, help='1D scanning z (default=False)') +parser.add_argument('-g', '--group', default=False, action='store_const', const=True, help='-g enables grouping') +parser.add_argument('--fill', default=False, action='store_const', const=True, help='Fill cavities after Grouping (only for 2D or 3D scanning) EXPERIMENTAL') + +parser.add_argument('--odx', default=False, help='OpenDX outfile') +parser.add_argument('--opdb', default=False, help='PDB outfile') + +parser.add_argument('--smd', default='0.00', help='Solvent minimum density (0.00 Atoms / A^3)') +parser.add_argument('--pmd', default='0.005', help='Protein minimum density (0.005 Atoms / A^3)') +parser.add_argument('--md', default='0', help='Minimum cavity diameter [A] (default=0)') +parser.add_argument('--overlap', default=False, action='store_const', const=True, help='Store ISV overlapping with protein voxels, FALSE per default') + + +parser.add_argument('--mdk', default=False,action='store_const',const=True, help='Minimum cavity diameter [A] (defaul=False) groups ISV > "--md [A] "and ISV < "--md [A]" separately --g MUST be set, currently only support with --sm 2D EXPERIMENTAL!') + +parser.add_argument('--fn', default=False,choices=range(1,27), help='Apply a neighbor filter int (1-26) (default=false)') + +args = parser.parse_args() + +''' +print args.solvent #solvent dx file +print args.protein # protein dxfile +print args.pmd # protein min density +print args.smd #solvent min density +print args.sm # scann method {1D 2D 3D) +print args.sm_x # enable scanning along x axis (same for y and z) +print args.group # enable grouping +print args.odx # output dxfile +print args.opdb # output pdb file +print args.md # min diameter of the cavity +print args.fn # filter by neighbor +args.overlap # switches if ISV overlapping with protein voxels > proteinMinDens will be deleted during scanning +''' + +print "dxTuber v" + version + +### +# Parameter check block :] +### + + + +###### +# checking if args.sm was given correct +###### +if not re.match('1D', args.sm): + if (args.sm_x or args.sm_y or args.sm_x): +# running 'dxtuber_cmd -h' : + sysRun = subprocess.Popen([sys.argv[0], '-h']) + sysRun.communicate() + print "Please use --sm 1D if --sm_x or --sm_y or sm_z are enabled" + sys.exit() + +if re.match('1D', args.sm): + if not (args.sm_x or args.sm_y or args.sm_z): + sysRun = subprocess.Popen([sys.argv[0], '-h']) + sysRun.communicate() + print "\nChoose at least one axis for 1D scanning:\n--sm 1D --sm_x and/or --sm_y and/or --sm_z" + + sys.exit() + +#### +# Cheking input files ... +#### + +if (args.protein or args.solvent) and args.ipdb: + sysRun = subprocess.Popen([sys.argv[0], '-h']) + sysRun.communicate() + print "\n Only Protein and Solvent OpenDX files (a) or cavities as pdb file (b) are supported as input file(s)\nCombinitions of (a) and (b) are not supported" + sys.exit() + +if not ((args.protein and args.solvent) or args.ipdb): + sysRun = subprocess.Popen([sys.argv[0], '-h']) + sysRun.communicate() + + print "\nInput Error:\nEither protein and solvent OpenDX files or a single pdb file are supported as input\n" + + sys.exit() + +#### +# fill + ipdb are not supported sind no solvent & protein densities are available in pdb format ... and yes they are nessacary for "fill" option +#### + +if args.ipdb and args.fill: + sysRun = subprocess.Popen([sys.argv[0], '-h']) + sysRun.communicate() + + print "\nFill option is not supported for pdb files as input\n" + sys.exit() + +###### +# Testing if any output will be created +###### + +if not (args.odx) and not (args.opdb): + sysRun = subprocess.Popen([sys.argv[0], '-h']) + sysRun.communicate() + + print "At least one output should be generated: \n--odx OUTPUT.DX or --opdb OUTPUT.PDB" + + sys.exit() + + +#### +# Testing args.fill +#### +if args.fill and args.sm == '1D' or (args.fill and not args.group): + sysRun = subprocess.Popen([sys.argv[0], '-h']) + sysRun.communicate() + print "Fill after grouping routine works only with 2D or 3D scanning, grouping is mandatory" + + sys.exit() + +if args.fill and args.mdk: + sysRun = subprocess.Popen([sys.argv[0], '-h']) + sysRun.communicate() + print "Fill after grouping in combinition with --mdk is not supported, ... not yet" + + sys.exit() +################## +# Oke parameters are fine, lets start :> +################## + + +if args.protein and args.solvent: + #### + #Read protein and water densities. + #### + + readDx = ReadOpenDx (args.protein) + readDx.importDx() + proteinDxObject = readDx.getDxObject() + readDx.setDxfile(args.solvent) + readDx.importDx() + solventDxObject = readDx.getDxObject() + + + ##################### + #Scanning + ##################### +# array of array of densities + scannedArray = False + + scanType = args.sm + scanDx = PlaneSearch (dxObjectWater = solventDxObject, dxObjectProtein = proteinDxObject, minDiameter=args.md,minProteinDens =args.pmd,minSolventDens = args.smd, mdk = args.mdk, storeOverLappingISV = args.overlap) + + + if args.sm == '1D': + scanMethod = 1 + if args.sm_x: + scanDx.oneDimension ('x') + scanType += " x" + if args.sm_y: + scanDx.oneDimension ('y') + scanType += " y" + if args.sm_z: + scanDx.oneDimension ('z') + scanType += " z" + scanned = scanDx.getDxObject() + if args.sm == '2D': + scanMethod = 2 + if args.mdk: + scanDx.twoDimensionMinDiaArray() + scannedArray = scanDx.getDxObjectArray() +# print str(len(scannedArray)) + else: + scanDx.twoDimension() + + scanned = scanDx.getDxObject() + if args.sm == '3D': + scanMethod = 3 + scanDx.threeDimension() + scanned = scanDx.getDxObject() + if args.mdk: + scannedArray = scanDx.getDxObjectArray() + +if args.ipdb: + readPDB = ReadPDB(args.ipdb) + readPDB.readPDB() + + scanned = readPDB.getDxObject() + scanMethod = readPDB.getTubeDxObject().getScanMethod() + scanType = readPDB.getTubeDxObject().getScanned() + args.pmd =readPDB.getTubeDxObject().getProtThreshold() + args.smd =readPDB.getTubeDxObject().getSolventThreshold() + filterName = readPDB.getTubeDxObject().getFilterApplied() + groupes = readPDB.getTubeDxObject().getGroupes() + grouped = readPDB.getTubeDxObject().getGrouped() + voxelGroupGrid = readPDB.getTubeDxObject().getVoxelGroupGrid() + args.protein = readPDB.getTubeDxObject().getProtFile() + args.solvent = readPDB.getTubeDxObject().getSolvFile() +################ +# Filter +################ +if args.protein and args.solvent: + filterName = 'None' + +if args.fn: + filter = FilterDx (scanned, args.md) + filter.filterTunnelNeighbour(neighbourMin=args.fn) + scanned = filter.getDxObject() + filterName += 'Neighbor ' + args.fn + +####### +# Grouping +####### +if args.protein and args.solvent: + grouped = 'No' + +if args.group: + if args.mdk: +## +# grouping isv > args.md +### + groupes = False + groupDx= GroupTubes(scannedArray[0]) + groupDx.findGroups() + groupesArray = [] + voxelGroupGridArray = [] + groupesArray.append(groupDx.getTubeGroup()) # 0 + voxelGroupGridArray.append(groupDx.getGroups()) #0 +# voxelGroupGridArray.extend(groupDx.getGroups()) #0 + +# they only contain the filtered ISV's + voxelGroupGrid = groupDx.getTubeGroup() + grouped = groupDx.getGroups() + groupDx= "" +### +# grouping isv < args.md +### + groupDx= GroupTubes (scannedArray[1]) + groupDx.findGroups() + groupesArray.append(groupDx.getTubeGroup()) # 1 + voxelGroupGridArray.append(groupDx.getGroups()) #1 #<<<<<----- these items should be sorted .... into one array ? [x][y][z] = density ? + #voxelGroupGridArray.extend(groupDx.getGroups()) #1 + grouped = 'Yes' + else: + groupesArray = False + voxelGroupGridArray = False + scannedArray = False + groupDx= GroupTubes (scanned) + groupDx.findGroups() + groupes = groupDx.getTubeGroup() + voxelGroupGrid = groupDx.getGroups() + grouped = 'Yes' + +### +v0.12 fill after grouping should be testet if it works in combinition with mdk *scared* + + if args.fill: + print "fill" + fill = FillTube(dxObjectWater = solventDxObject, + dxObjectProtein = proteinDxObject, + dxObjectTube = scanned, + scanMethod = scanMethod, + minProteinDensity = args.pmd, + minSolvDensity = args.smd, + groups = groupes, + voxelGroupGrid = voxelGroupGrid + ) + fill.fillGroups() + voxelGroupGrid=fill.getVoxelGroupGrid() + groupes=fill.getGroups() + scanned = fill.getNewTubeObject() +### filltube needs to be rewritten really ugly code inside .... just to be sure nothing strange happens group again ... + groupDx= GroupTubes (scanned) + groupDx.findGroups() + groupes = groupDx.getTubeGroup() + voxelGroupGrid = groupDx.getGroups() +### / fill after grouping + + result = TubeDx ( scanned, + filename="None", + scanMethod = scanMethod, + scanned = scanType, + solventThreshold = args.smd , + protThreshold = args.pmd, + minDiameter = args.md, + protFile = args.protein, + solvFile = args.solvent, + version = version, + filterApplied = filterName, + groupes = groupes, + VoxelGroupGrid = voxelGroupGrid, + grouped = grouped, + groupesArray = groupesArray, + VoxelGroupGridArray = voxelGroupGridArray, + dxObjectArray = scannedArray + ) + + + +### else ... no grouping flag ... +else: + if args.protein and args.solvent: + result = TubeDx ( scanned, + filename="None", + scanMethod = scanMethod, + scanned = scanType, + solventThreshold = args.smd , + protThreshold = args.pmd, + minDiameter = args.md, + protFile = args.protein, + solvFile = args.solvent, + version = version, + filterApplied = filterName, + # groupes = groupes, + # VoxelGroupGrid = voxelGroupGrid, + grouped = grouped + ) + if args.ipdb: + result = readPDB.getTubeDxObject() + + + +########### +# Output +########### + + +if args.odx: + saveToDx = WriteDx (scanned) + saveToDx.write(args.odx) + +if args.opdb: + saveToPDB = WritePDB (result) + if args.group: + if args.mdk: +# print result.getGroupesArray() + saveToPDB.writeMultiGroups(result.getGroupesArray(),args.opdb) + else: + saveToPDB.writeGroups(result.getGroupes(),args.opdb) + + else: + saveToPDB.write(args.opdb) + +''' +storePDB = WritePDB (result) +#saveToDx = WriteDx (scanned) +outfile = sys.argv[1][:-3] +outfile += ".pdb" +#outfile += "_out.dx" +storePDB.writeGroups(result.getGroupes(),outfile) +#saveToDx.write(outfile) +''' + + + + diff --git a/dxTuber_gtk.py b/dxTuber_gtk.py new file mode 100755 index 0000000..235eeeb --- /dev/null +++ b/dxTuber_gtk.py @@ -0,0 +1,474 @@ +#!/usr/bin/env python +# +# This file is part of the dxTuber package +# dxTuber v0.27 7. Jan 2013 +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +# +#################################################################################### +### This tool s written to detect cavities based on protein and solvent movement ### +### ### +### ### +### 1. flood your system with water (by using GROMACS) ### +### 2. Open your molecule in VMD ### +### 3. save 2 OpenDX (only solvent and only protein) files via VMD -> VolMap ### +### 4. start this script ### +### ### +### Contact: raunest@bit.uni-bonn.de or m.raunest@gmx.de Martin Raunest ### +### ### +#################################################################################### + + +import sys +import os +import re +####################### +# only within gnome nessasary KDE ? +####################### +import gnome.ui +gnome.init("dxTuber", "version") +####################### +import pygtk +try: + import gtk + import gtk.glade +except: + print "You need to install pyGTK or GTKv2 ", + print "or set your PYTHONPATH correctly." + print "try: export PYTHONPATH=", + print "/usr/local/lib/python2.X/site-packages/" + sys.exit(1) + +from lib.readopendx import * +from lib.planesearch import * +from lib.writepdb import * +from lib.writedx import * +from lib.grouptubes import * +from lib.filterdx import * +from lib.analysetube import * +from lib.filltube import * +from lib.tubedx import * +from lib.readpdb import * + +class appgui: + __dxTubeArray = [] + __proteinSelected = 0 + __waterSelected = 0 + __fileSelected = 0 + __scanMethod = 2 + __analyseAxis=0 + __analyseMethod = 0 + def toggleProteinCB(self,widget): + list = widget.get_group() + list.reverse() + for button in range (len(list)): + if (list[button].get_active()): + self.__proteinSelected = button + + def toggleWaterCB(self,widget): + list = widget.get_group() + list.reverse() + for button in range (len(list)): + if (list[button].get_active()): + self.__waterSelected = button + + def toggleFileCB(self,widget): + list = widget.get_group() + list.reverse() + for button in range (len(list)): + if (list[button].get_active()): + self.__fileSelected = button +### Scan Methods +#1: 1D +#2: 2D +#3: 3D +### + def toggleScan1D(self,widget): + xButton = self.wTree.get_widget("search1DXCheckBox") + yButton = self.wTree.get_widget("search1DYCheckBox") + zButton = self.wTree.get_widget("search1DZCheckBox") + + xButton.set_sensitive(1) + yButton.set_sensitive(1) + zButton.set_sensitive(1) + self.__scanMethod = 1 + def toggleScan2D(self,widget): + xButton = self.wTree.get_widget("search1DXCheckBox") + yButton = self.wTree.get_widget("search1DYCheckBox") + zButton = self.wTree.get_widget("search1DZCheckBox") + + xButton.set_sensitive(0) + yButton.set_sensitive(0) + zButton.set_sensitive(0) + self.__scanMethod = 2 + def toggleScan3D(self,widget): + xButton = self.wTree.get_widget("search1DXCheckBox") + yButton = self.wTree.get_widget("search1DYCheckBox") + zButton = self.wTree.get_widget("search1DZCheckBox") + + xButton.set_sensitive(0) + yButton.set_sensitive(0) + zButton.set_sensitive(0) + self.__scanMethod = 3 + + def toggleAnalyseAxisCB(self,widget): + list = widget.get_group() + list.reverse() + for axis in range (len(list)): + if (list[axis].get_active()): + self.__analyseAxis = axis + def toggleAnalyseMethodCB(self,widget): + list = widget.get_group() + list.reverse() + for method in range (len(list)): + if (list[method].get_active()): + self.__analyseMethod = method + + def __readSettings(self): +# self.__version = '0.09' +# self.__vDate = '17.Mar.2011' + dxTuberConfFile = os.environ['HOME']+"/.dxTuber/dxTuber.conf" + conf = open (dxTuberConfFile,'r') + confLines = conf.readlines() + for line in confLines: + if re.match ('^#',line): + continue + if re.match('INSTALL_DIR',line): + tmpArray = re.split('\s',line,1) + tmp = re.sub('\"','',tmpArray[1]) + tmp = re.sub('\n','',tmp) + self.instPath = re.sub('^\s','',tmp) + if re.match('VERSION',line): +# print line + tmpArray = re.split('\s',line,1) + tmp = re.sub('\"','',tmpArray[1]) + tmp = re.sub('\n','',tmp) + self.__version = re.sub('^\s','',tmp) + print "dxTuber v" + self.__version + + def __init__(self): + self.__readSettings() + glade = self.instPath+"/dxTuber.libglade" + self.wTree=gtk.glade.XML (fname=str(glade), domain="app1") + self.__tableWidget = self.wTree.get_widget("newTable") + scanButton = self.wTree.get_widget("scanButton") + scanButton.connect("clicked", self.scanTube) + + dic = { "py_open_file" : self.openDx, + "OPEN_PDB" : self.openPDB, + "FILE_CHOOSER_ACTION_SAVE": self.savePDB, + "SAVE_TO_DX" : self.saveDx, + "GROUPE_DX" : self.groupTube, + "SCAN_1D_CLICKED" : self.toggleScan1D, + "SCAN_2D_CLICKED" : self.toggleScan2D, + "SCAN_3D_CLICKED" : self.toggleScan3D, + "ANALYSE_TUBE_BUTTON_CLICKED" :self.startAnalysePopUp, + "FILTER_TUBE_BUTTON_CLICKED": self.filterTube, + } + self.wTree.signal_autoconnect(dic) + def filterTube(self,widget): + filterNeighbourRadio = self.wTree.get_widget("filterNeighbourRadio") + filterDistanceSpinBox = self.wTree.get_widget("filterDistanceSpinBox") + filter = FilterDx (self.__dxTubeArray[self.__fileSelected].getDxObject(), filterDistanceSpinBox.get_value_as_int()) + if filterNeighbourRadio.get_active(): + filterSpinBox = self.wTree.get_widget("filterNeighbourSpinbox") + filterValue = filterSpinBox.get_value_as_int() + filter.filterTunnelNeighbour(neighbourMin=filterValue) #parameter of neighbours should be set in future + filterName = 'Neighbor ' + str(filterSpinBox.get_value_as_int()) + else: + filter.filterDistance() + filterName = 'Distance '+ str(filterDistanceSpinBox.get_value_as_int()) + + self.__dxTubeArray.append (TubeDx(filter.getDxObject(),filename="None", + scanMethod = self.__dxTubeArray[self.__fileSelected].getScanMethod(), + scanned = self.__dxTubeArray[self.__fileSelected].getScanned(), + solventThreshold = self.__dxTubeArray[self.__fileSelected].getSolventThreshold(), + protThreshold = self.__dxTubeArray[self.__fileSelected].getProtThreshold(), + minDiameter = self.__dxTubeArray[self.__fileSelected].getMinDiameter(), + protFile = self.__dxTubeArray[self.__fileSelected].getProtFile(), + solvFile = self.__dxTubeArray[self.__fileSelected].getSolvFile(), + version = self.__version, + filterApplied = self.__dxTubeArray[self.__fileSelected].getFilterApplied() + ) + ) + self.__dxTubeArray[(len(self.__dxTubeArray)-1)].setFilterApplied(filterName) + print "Done" + self.gtkDrawFileTable() + def scanTube (self,widget = None): + proteinDx = self.__dxTubeArray[self.__proteinSelected].getDxObject() + waterDx = self.__dxTubeArray[self.__waterSelected].getDxObject() + protThreshEntry = self.wTree.get_widget("proteinThresholdEntry") +###### implement solvent densities and min diameter: + solventThreshEntry = self.wTree.get_widget("solventThresholdEntry") + minDiameterEntry = self.wTree.get_widget("minDiameterEntry") + scanDx = PlaneSearch (dxObjectWater = waterDx, dxObjectProtein = proteinDx, + minDiameter = minDiameterEntry.get_text(), + minProteinDens = protThreshEntry.get_text(), + minSolventDens = solventThreshEntry.get_text() + ) + oneDimRadio = self.wTree.get_widget("search1DRadio") + if self.__scanMethod == 1: + scanType = "1D" + xButton = self.wTree.get_widget("search1DXCheckBox") + yButton = self.wTree.get_widget("search1DYCheckBox") + zButton = self.wTree.get_widget("search1DZCheckBox") + if xButton.get_active(): + scanDx.oneDimension ('x') + scanType += " x" + if yButton.get_active(): + scanDx.oneDimension ('y') + scanType += " y" + if zButton.get_active(): + scanDx.oneDimension ('z') + scanType += " z" + if self.__scanMethod == 2: + scanDx.twoDimension() + scanType = "2D" + if self.__scanMethod == 3: + scanDx.threeDimension() + scanType = "3D" + filename = 'None' + self.__dxTubeArray.append (TubeDx(scanDx.getDxObject(),filename,scanMethod = self.__scanMethod,scanned =scanType,waterObject = self.__dxTubeArray[self.__waterSelected].getDxObject(), proteinObject=self.__dxTubeArray[self.__proteinSelected].getDxObject(),minDiameter = minDiameterEntry.get_text() ,protThreshold = protThreshEntry.get_text(), solventThreshold = solventThreshEntry.get_text(), protFile = proteinDx.getFilename(), solvFile = waterDx.getFilename(), version = self.__version)) + self.gtkDrawFileTable() + + def startAnalysePopUp(self,widget): + glade = self.instPath+"/analyseDialog.libglade" + print glade + self.analyseDialog=gtk.glade.XML (fname=str(glade), domain="analyseDialog") + analysePopup = self.analyseDialog.get_widget("analysePopup") + tablePopup = self.analyseDialog.get_widget("table1") + + tubeObject = self.__dxTubeArray[self.__fileSelected] + self.groupeCombobox = gtk.combo_box_new_text() + for tube in range (len (tubeObject.getGroupes() )): + self.groupeCombobox.append_text(str(tube)) + tablePopup.attach(self.groupeCombobox,1,2,1,2,yoptions=0) + + analysePopup.show() + self.groupeCombobox.show() + dic = { "ANALYSE_TUBE_BUTTON_CLICKED" : self.startAnalyse, + "ANALYSE_X" : self.toggleAnalyseAxisCB, + "ANALYSE_Y" : self.toggleAnalyseAxisCB, + "ANALYSE_Z" : self.toggleAnalyseAxisCB, + "ANALYSE_METHOD_0": self.toggleAnalyseMethodCB, + "ANALYSE_METHOD_1": self.toggleAnalyseMethodCB + } + self.analyseDialog.signal_autoconnect(dic) + + def startAnalyse (self,widget=None): + tubeCounter = int( self.groupeCombobox.get_active() ) +#voxelGroupGrid +#[x][y][z] = group + analyse = AnalyseTube (self.__dxTubeArray[self.__fileSelected].getVoxelGroupGrid(),tubeCounter) + chooser = gtk.FileChooserDialog(title=None,action=gtk.FILE_CHOOSER_ACTION_SAVE,buttons=(gtk.STOCK_CANCEL,gtk.RESPONSE_CANCEL,gtk.STOCK_SAVE,gtk.RESPONSE_OK)) + chooser.run() + filename = chooser.get_filename() +# method 0 is area + if self.__analyseMethod == 0: + if self.__analyseAxis == 0: + analyse.analyseTubeArea('x') + analyse.writeStatistics(filename) + if self.__analyseAxis == 1: + analyse.analyseTubeArea('y') + analyse.writeStatistics(filename) + if self.__analyseAxis == 2: + analyse.analyseTubeArea('z') + analyse.writeStatistics(filename) +#method 1 is max diameter + if self.__analyseMethod == 1: + if self.__analyseAxis == 0: + analyse.analyseTubeDiameterMax('x') + analyse.writeStatistics(filename) + if self.__analyseAxis == 1: + analyse.analyseTubeDiameterMax('y') + analyse.writeStatistics(filename) + if self.__analyseAxis == 2: + analyse.analyseTubeDiameterMax('z') + analyse.writeStatistics(filename) + chooser.destroy() + self.analyseDialog.get_widget("analysePopup").destroy() + def groupTube (self,widget = None): + dxTube = self.__dxTubeArray[self.__fileSelected] + tubes = GroupTubes (dxTube.getDxObject()) + tubes.findGroups() + dxTube.setGroupes (tubes.getTubeGroup()) + dxTube.setVoxelGroupGrid (tubes.getGroups()) + fillTubeButton = self.wTree.get_widget("fillTubeCheckButton") + if fillTubeButton.get_active(): + protThreshEntry = self.wTree.get_widget("proteinThresholdEntry") + protThresh = protThreshEntry.get_text() + solvThreshEntry = self.wTree.get_widget("solventThresholdEntry") + solvThresh = solvThreshEntry.get_text() + fill = FillTube(dxObjectWater = dxTube.getWaterObject(), + dxObjectProtein = dxTube.getProteinObject(), + dxObjectTube = dxTube.getDxObject(), + scanMethod = dxTube.getScanMethod(), + minProteinDensity = protThresh, + minSolvDensity = solvThresh, + groups = dxTube.getGroupes(), + voxelGroupGrid = dxTube.getVoxelGroupGrid() + ) + fill.fillGroups() + dxTube.setVoxelGroupGrid(fill.getVoxelGroupGrid()) + dxTube.setGroupes(fill.getGroups()) + dxTube.setDxObject(fill.getNewTubeObject()) + print "Finished cavity filling" + dxTube.setGrouped("Yes") + self.gtkDrawFileTable() + def openDx (self,file): + chooser = gtk.FileChooserDialog(title=None,action=gtk.FILE_CHOOSER_ACTION_OPEN, + buttons=(gtk.STOCK_CANCEL,gtk.RESPONSE_CANCEL,gtk.STOCK_OPEN,gtk.RESPONSE_OK)) + chooser.run() + filename = chooser.get_filename() + chooser.destroy() + readDx = ReadOpenDx (filename) + readDx.importDx() + self.__dxTubeArray.append (TubeDx(readDx,filename)) + self.gtkDrawFileTable() + def openPDB (self, widget = None): + chooser = gtk.FileChooserDialog(title=None,action=gtk.FILE_CHOOSER_ACTION_OPEN, + buttons=(gtk.STOCK_CANCEL,gtk.RESPONSE_CANCEL,gtk.STOCK_OPEN,gtk.RESPONSE_OK)) + chooser.run() + filename = chooser.get_filename() + chooser.destroy() + readPDB = ReadPDB ( filename) + if readPDB.readPDB() == 0: + self.__dxTubeArray.append (readPDB.getTubeDxObject()) + self.gtkDrawFileTable() + def saveDx(self,widget=None): + dxObject = self.__dxTubeArray[self.__fileSelected].getDxObject() + chooser = gtk.FileChooserDialog(title=None,action=gtk.FILE_CHOOSER_ACTION_SAVE,buttons=(gtk.STOCK_CANCEL,gtk.RESPONSE_CANCEL,gtk.STOCK_SAVE,gtk.RESPONSE_OK)) + chooser.run() + filename = chooser.get_filename() + chooser.destroy() + saveToDx = WriteDx (dxObject) + saveToDx.write(filename) + if ( re.match("^None",self.__dxTubeArray[self.__fileSelected].getFilename())): + self.__dxTubeArray[self.__fileSelected].setFilename(filename) + self.gtkDrawFileTable() + def savePDB(self,file): + chooser = gtk.FileChooserDialog(title=None,action=gtk.FILE_CHOOSER_ACTION_SAVE,buttons=(gtk.STOCK_CANCEL,gtk.RESPONSE_CANCEL,gtk.STOCK_SAVE,gtk.RESPONSE_OK)) + chooser.run() + filename = chooser.get_filename() + temp = self.__dxTubeArray[self.__fileSelected] + savePDB = WritePDB (temp) + if self.__dxTubeArray[self.__fileSelected].getGrouped() == 'Yes': + savePDB.writeGroups(self.__dxTubeArray[self.__fileSelected].getGroupes(),filename) + else: + savePDB.write(filename) + chooser.destroy() + + def gtkDrawFileTable (self): +#this will also destroy the notebooklabel ! + self.__tableWidget.destroy() + noteWidget = self.wTree.get_widget("notebook1") + newTable = gtk.Table (rows = (len (self.__dxTubeArray)+1), columns= 9) + + noteLabel = gtk.Label (str="Files") + fileLabel = gtk.Label (str="Files") + setAsProtLabel = gtk.Label (str="Set as protein") + setAsWaterLabel = gtk.Label (str="Set as water") + setAsScannedLabel= gtk.Label (str="Scanned") + groupedLabel = gtk.Label (str="Grouped") + setCalcAvLabel = gtk.Label (str="Calc averages") + minDiameterLabel = gtk.Label (str="Min diameter") + minProtThresholdLabel = gtk.Label (str="Min protein density") + minSolventLabel = gtk.Label (str="Min solvent density") + + + noteWidget.insert_page(newTable, tab_label = noteLabel, position =0) + newTable.show () + noteLabel.show() + newTable.attach (fileLabel,0,1,0,1,yoptions=0) + newTable.attach (setAsProtLabel,1,2,0,1,yoptions=0) + newTable.attach (setAsWaterLabel,2,3,0,1,yoptions=0) + newTable.attach (setAsScannedLabel,3,4,0,1,yoptions=0) + newTable.attach (groupedLabel,4,5,0,1,yoptions=0) + newTable.attach (minDiameterLabel,5,6,0,1,yoptions=0) + newTable.attach (minProtThresholdLabel,6,7,0,1,yoptions=0) + newTable.attach (minSolventLabel,7,8,0,1,yoptions=0) +# newTable.attach (setCalcAvLabel,5,6,0,1,yoptions=0) + + setAsScannedLabel.show() + fileLabel.show() + setAsProtLabel.show() + setAsWaterLabel.show() + groupedLabel.show() + minDiameterLabel.show() + minProtThresholdLabel.show() + minSolventLabel.show() +# setCalcAvLabel.show() + selectProtein = self.__proteinSelected +# print selectProtein + for lineCounter in range (len (self.__dxTubeArray)): + if (lineCounter == 0): + self.__firstFileRadio = gtk.RadioButton (label=self.__dxTubeArray[lineCounter].getFilename(), use_underline=False) + newTable.attach (self.__firstFileRadio,0,1,lineCounter+1,lineCounter+2,yoptions=0) + self.__firstFileRadio.show() + self.__firstFileRadio.connect("toggled", self.toggleFileCB) + self.toggleFileCB (self.__firstFileRadio) + else: + fileRadio = gtk.RadioButton (group = self.__firstFileRadio, label=self.__dxTubeArray[lineCounter].getFilename(), use_underline=False) + newTable.attach (fileRadio,0,1,lineCounter+1,lineCounter+2,yoptions=0) + fileRadio.show() + fileRadio.connect("toggled", self.toggleFileCB) + if (lineCounter == 0): + self.__firstProtRadioButton = gtk.RadioButton () + newTable.attach (self.__firstProtRadioButton,1,2,lineCounter+1,lineCounter+2,yoptions=0) + self.__firstProtRadioButton.show() + self.__firstProtRadioButton.connect("toggled", self.toggleProteinCB) + self.toggleProteinCB(self.__firstProtRadioButton) + else: + setAsProtRadioButton = gtk.RadioButton (group = self.__firstProtRadioButton) + newTable.attach (setAsProtRadioButton,1,2,lineCounter+1,lineCounter+2,yoptions=0) + setAsProtRadioButton.show() + setAsProtRadioButton.connect("toggled", self.toggleProteinCB) + + if (lineCounter == 0): + self.__firstWaterRadioButton = gtk.RadioButton () + newTable.attach (self.__firstWaterRadioButton,2,3,lineCounter+1,lineCounter+2,yoptions=0) + self.__firstWaterRadioButton.show() + self.__firstWaterRadioButton.connect("toggled", self.toggleWaterCB) + self.toggleWaterCB(self.__firstWaterRadioButton) + else: + setAsWaterRadioButton = gtk.RadioButton (group = self.__firstWaterRadioButton) + newTable.attach (setAsWaterRadioButton,2,3,lineCounter+1,lineCounter+2,yoptions=0) + setAsWaterRadioButton.show() + setAsWaterRadioButton.connect("toggled", self.toggleWaterCB) + + groupedLabel = gtk.Label (str=self.__dxTubeArray[lineCounter].getGrouped()) + setAsScannedLabel= gtk.Label (str=self.__dxTubeArray[lineCounter].getScanned()) + minDiameterValue = gtk.Label (str = self.__dxTubeArray[lineCounter].getMinDiameter()) + minProtThresholdValue = gtk.Label (str = self.__dxTubeArray[lineCounter].getProtThreshold()) + minSolventValue = gtk.Label (str = self.__dxTubeArray[lineCounter].getSolventThreshold()) + + newTable.attach (setAsScannedLabel,3,4,lineCounter+1,lineCounter+2,yoptions=0) + newTable.attach (groupedLabel,4,5,lineCounter+1,lineCounter+2,yoptions=0) + newTable.attach (minDiameterValue,5,6,lineCounter+1,lineCounter+2,yoptions=0) + newTable.attach (minProtThresholdValue,6,7,lineCounter+1,lineCounter+2,yoptions=0) + newTable.attach (minSolventValue,7,8,lineCounter+1,lineCounter+2,yoptions=0) + + setAsScannedLabel.show() + groupedLabel.show() + minDiameterValue.show() + minProtThresholdValue.show() + minSolventValue.show() + + self.__tableWidget = newTable + noteWidget.set_current_page(0) +app=appgui() +gtk.main() + diff --git a/feature_request b/feature_request new file mode 100644 index 0000000..a3a5cc0 --- /dev/null +++ b/feature_request @@ -0,0 +1,4 @@ +dichten in tuber wieder ruawerfen aus der gui (unload) +known issues (bugs die bekannt sind) +spalten in der gui anpassen +springede radio buttons diff --git a/filterVoxel.py b/filterVoxel.py new file mode 100755 index 0000000..e95cbfa --- /dev/null +++ b/filterVoxel.py @@ -0,0 +1,165 @@ +#!/usr/bin/env python +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +### +# +# v0.2 +# +yay lets start posst proccessing +# + +import sys +import os +from lib.readopendx import * +from lib.writedx import * +from lib.writepdb import * +from lib.filterdx import * +from lib.readpdb import * +from lib.tubedx import * + +######### +# Get dxTuber settings +######### + + + +dxTuberConfFile = os.environ['HOME']+"/.dxTuber/dxTuber.conf" +conf = open (dxTuberConfFile,'r') +confLines = conf.readlines() +for line in confLines: + if re.match ('^#',line): + continue + if re.match('VERSION',line): + tmpArray = re.split('\s',line,1) + tmp = re.sub('\"','',tmpArray[1]) + tmp = re.sub('\n','',tmp) + version = re.sub('^\s','',tmp) + +#### +# Deal with arguments +#### + +###################################################################### +# If argparse is not installed global, get a static lib # +# and import it manually (static) # +# # +# http://packages.ubuntu.com/lucid/python-argparse # +# http://code.google.com/p/argparse/ # +# # +# maybe this helps for implementation ... # +# # +# sys.path.append('path_to_argparse_lib') +# from argparse import * +# # +###################################################################### + + +import argparse + +parser = argparse.ArgumentParser(description='FilterVoxel was written to post proccess results of dxTuber:\nMin and max defines coordinates of an rectangle',epilog='dxTuber v' + version) + +### +# for static libs use: +# parser = ArgumentParser(description='Command line version of dxTuber',epilog='dxTuber v' + version) +### +parser.add_argument('--idx',default=False,help='Input OpenDX file') +parser.add_argument('--ipdb',default=False,help='Input PDB file') +parser.add_argument('--odx',default=False,help='Output OpenDX file') +parser.add_argument('--opdb',default=False,help='Output PDB file') +parser.add_argument('-g', '--group', default=False, action='store_const', const=True, help='-g enables grouping') + +parser.add_argument('--min',default=False,help='Array containing minimum "x y z" coordinates') +parser.add_argument('--max',default=False,help='Array containing maximum "x y z" coordinates') + + + +args = parser.parse_args() + +print args.min + +if (not args.idx and not args.ipdb) or (args.idx and args.ipdb) : + print "Define ONE input file either ipdb OR idx" + sys.exit() + +if args.group and args.idx: + print "Goupes can not be stored in OpenDX file format" + sys.exit() + +minArray = args.min.split() +maxArray = args.max.split() + +if args.idx: + readDx = ReadOpenDx (args.idx) + readDx.importDx() + dxObject = readDx.getDxObject() +if args.pdb: + readPDB = ReadPDB (args.pdb) + readPDB.readPDB() + dxObject = readPDB.getDxObject() + +filterThis = FilterDx (dxObject) +filterThis.deleteRectangleInside(minArray,maxArray) +dxObject = filterThis.getDxObject() + + +if args.group: + groupDx= GroupTubes (dxObject) + groupDx.findGroups() + groupes = groupDx.getTubeGroup() + voxelGroupeGrid = groupDx.getGroups() +# grouped = 'Yes' + +if args.odx: + saveToDx = WriteDx (dxObject) + saveToDx.write(args.odx) + +if args.opdb: + saveToPDB = WritePDB (result) + if args.group: + if args.ipdb: + tubeDxObject = readDx.getTubeDxObject() + tubeDxObject.setDxObject(dxObject) + if tubeDxObject.getFilterApplied() == 'None': + tubeDxObject.setFilterApplied("delRectangleInside") + else: + tubeDxObject.setFilterApplied(tubeDxObject.getFilterApplied() + " delRectangleInside") + tubeDxObject.setVoxelGroupGrid(voxelGroupeGrid) + tubeDxObject.setGroupes(groupes) + tubeDxObject.setGrouped("Yes") + tubeDxObject.setFilename(args.ipdb) + + if args.idx: + tubeDxObject = TubeDx ( dxObject, + filename=args.idx, +# scanMethod = "n/a", +# scanned = "n/a", +# solventThreshold = "n/a" , +# protThreshold = "n/a", #better keep theire default values ... +# minDiameter = args.md, +# protFile = args.protein, +# solvFile = args.solvent, + version = version, + filterApplied = "delRectangleInside", + groupes = groupes, + VoxelGroupGrid = voxelGroupeGrid, + grouped = 'Yes' + ) + + + else: + saveToPDB.write(args.opdb) + + diff --git a/gpl-2.0.txt b/gpl-2.0.txt new file mode 100644 index 0000000..d159169 --- /dev/null +++ b/gpl-2.0.txt @@ -0,0 +1,339 @@ + GNU GENERAL PUBLIC LICENSE + Version 2, June 1991 + + Copyright (C) 1989, 1991 Free Software Foundation, Inc., + 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA + Everyone is permitted to copy and distribute verbatim copies + of this license document, but changing it is not allowed. + + Preamble + + The licenses for most software are designed to take away your +freedom to share and change it. By contrast, the GNU General Public +License is intended to guarantee your freedom to share and change free +software--to make sure the software is free for all its users. This +General Public License applies to most of the Free Software +Foundation's software and to any other program whose authors commit to +using it. (Some other Free Software Foundation software is covered by +the GNU Lesser General Public License instead.) You can apply it to +your programs, too. + + When we speak of free software, we are referring to freedom, not +price. Our General Public Licenses are designed to make sure that you +have the freedom to distribute copies of free software (and charge for +this service if you wish), that you receive source code or can get it +if you want it, that you can change the software or use pieces of it +in new free programs; and that you know you can do these things. + + To protect your rights, we need to make restrictions that forbid +anyone to deny you these rights or to ask you to surrender the rights. +These restrictions translate to certain responsibilities for you if you +distribute copies of the software, or if you modify it. + + For example, if you distribute copies of such a program, whether +gratis or for a fee, you must give the recipients all the rights that +you have. You must make sure that they, too, receive or can get the +source code. And you must show them these terms so they know their +rights. + + We protect your rights with two steps: (1) copyright the software, and +(2) offer you this license which gives you legal permission to copy, +distribute and/or modify the software. + + Also, for each author's protection and ours, we want to make certain +that everyone understands that there is no warranty for this free +software. If the software is modified by someone else and passed on, we +want its recipients to know that what they have is not the original, so +that any problems introduced by others will not reflect on the original +authors' reputations. + + Finally, any free program is threatened constantly by software +patents. We wish to avoid the danger that redistributors of a free +program will individually obtain patent licenses, in effect making the +program proprietary. To prevent this, we have made it clear that any +patent must be licensed for everyone's free use or not licensed at all. + + The precise terms and conditions for copying, distribution and +modification follow. + + GNU GENERAL PUBLIC LICENSE + TERMS AND CONDITIONS FOR COPYING, DISTRIBUTION AND MODIFICATION + + 0. This License applies to any program or other work which contains +a notice placed by the copyright holder saying it may be distributed +under the terms of this General Public License. The "Program", below, +refers to any such program or work, and a "work based on the Program" +means either the Program or any derivative work under copyright law: +that is to say, a work containing the Program or a portion of it, +either verbatim or with modifications and/or translated into another +language. (Hereinafter, translation is included without limitation in +the term "modification".) Each licensee is addressed as "you". + +Activities other than copying, distribution and modification are not +covered by this License; they are outside its scope. The act of +running the Program is not restricted, and the output from the Program +is covered only if its contents constitute a work based on the +Program (independent of having been made by running the Program). +Whether that is true depends on what the Program does. + + 1. You may copy and distribute verbatim copies of the Program's +source code as you receive it, in any medium, provided that you +conspicuously and appropriately publish on each copy an appropriate +copyright notice and disclaimer of warranty; keep intact all the +notices that refer to this License and to the absence of any warranty; +and give any other recipients of the Program a copy of this License +along with the Program. + +You may charge a fee for the physical act of transferring a copy, and +you may at your option offer warranty protection in exchange for a fee. + + 2. You may modify your copy or copies of the Program or any portion +of it, thus forming a work based on the Program, and copy and +distribute such modifications or work under the terms of Section 1 +above, provided that you also meet all of these conditions: + + a) You must cause the modified files to carry prominent notices + stating that you changed the files and the date of any change. + + b) You must cause any work that you distribute or publish, that in + whole or in part contains or is derived from the Program or any + part thereof, to be licensed as a whole at no charge to all third + parties under the terms of this License. + + c) If the modified program normally reads commands interactively + when run, you must cause it, when started running for such + interactive use in the most ordinary way, to print or display an + announcement including an appropriate copyright notice and a + notice that there is no warranty (or else, saying that you provide + a warranty) and that users may redistribute the program under + these conditions, and telling the user how to view a copy of this + License. (Exception: if the Program itself is interactive but + does not normally print such an announcement, your work based on + the Program is not required to print an announcement.) + +These requirements apply to the modified work as a whole. If +identifiable sections of that work are not derived from the Program, +and can be reasonably considered independent and separate works in +themselves, then this License, and its terms, do not apply to those +sections when you distribute them as separate works. But when you +distribute the same sections as part of a whole which is a work based +on the Program, the distribution of the whole must be on the terms of +this License, whose permissions for other licensees extend to the +entire whole, and thus to each and every part regardless of who wrote it. + +Thus, it is not the intent of this section to claim rights or contest +your rights to work written entirely by you; rather, the intent is to +exercise the right to control the distribution of derivative or +collective works based on the Program. + +In addition, mere aggregation of another work not based on the Program +with the Program (or with a work based on the Program) on a volume of +a storage or distribution medium does not bring the other work under +the scope of this License. + + 3. You may copy and distribute the Program (or a work based on it, +under Section 2) in object code or executable form under the terms of +Sections 1 and 2 above provided that you also do one of the following: + + a) Accompany it with the complete corresponding machine-readable + source code, which must be distributed under the terms of Sections + 1 and 2 above on a medium customarily used for software interchange; or, + + b) Accompany it with a written offer, valid for at least three + years, to give any third party, for a charge no more than your + cost of physically performing source distribution, a complete + machine-readable copy of the corresponding source code, to be + distributed under the terms of Sections 1 and 2 above on a medium + customarily used for software interchange; or, + + c) Accompany it with the information you received as to the offer + to distribute corresponding source code. (This alternative is + allowed only for noncommercial distribution and only if you + received the program in object code or executable form with such + an offer, in accord with Subsection b above.) + +The source code for a work means the preferred form of the work for +making modifications to it. For an executable work, complete source +code means all the source code for all modules it contains, plus any +associated interface definition files, plus the scripts used to +control compilation and installation of the executable. However, as a +special exception, the source code distributed need not include +anything that is normally distributed (in either source or binary +form) with the major components (compiler, kernel, and so on) of the +operating system on which the executable runs, unless that component +itself accompanies the executable. + +If distribution of executable or object code is made by offering +access to copy from a designated place, then offering equivalent +access to copy the source code from the same place counts as +distribution of the source code, even though third parties are not +compelled to copy the source along with the object code. + + 4. You may not copy, modify, sublicense, or distribute the Program +except as expressly provided under this License. Any attempt +otherwise to copy, modify, sublicense or distribute the Program is +void, and will automatically terminate your rights under this License. +However, parties who have received copies, or rights, from you under +this License will not have their licenses terminated so long as such +parties remain in full compliance. + + 5. You are not required to accept this License, since you have not +signed it. However, nothing else grants you permission to modify or +distribute the Program or its derivative works. These actions are +prohibited by law if you do not accept this License. Therefore, by +modifying or distributing the Program (or any work based on the +Program), you indicate your acceptance of this License to do so, and +all its terms and conditions for copying, distributing or modifying +the Program or works based on it. + + 6. Each time you redistribute the Program (or any work based on the +Program), the recipient automatically receives a license from the +original licensor to copy, distribute or modify the Program subject to +these terms and conditions. You may not impose any further +restrictions on the recipients' exercise of the rights granted herein. +You are not responsible for enforcing compliance by third parties to +this License. + + 7. If, as a consequence of a court judgment or allegation of patent +infringement or for any other reason (not limited to patent issues), +conditions are imposed on you (whether by court order, agreement or +otherwise) that contradict the conditions of this License, they do not +excuse you from the conditions of this License. If you cannot +distribute so as to satisfy simultaneously your obligations under this +License and any other pertinent obligations, then as a consequence you +may not distribute the Program at all. For example, if a patent +license would not permit royalty-free redistribution of the Program by +all those who receive copies directly or indirectly through you, then +the only way you could satisfy both it and this License would be to +refrain entirely from distribution of the Program. + +If any portion of this section is held invalid or unenforceable under +any particular circumstance, the balance of the section is intended to +apply and the section as a whole is intended to apply in other +circumstances. + +It is not the purpose of this section to induce you to infringe any +patents or other property right claims or to contest validity of any +such claims; this section has the sole purpose of protecting the +integrity of the free software distribution system, which is +implemented by public license practices. Many people have made +generous contributions to the wide range of software distributed +through that system in reliance on consistent application of that +system; it is up to the author/donor to decide if he or she is willing +to distribute software through any other system and a licensee cannot +impose that choice. + +This section is intended to make thoroughly clear what is believed to +be a consequence of the rest of this License. + + 8. If the distribution and/or use of the Program is restricted in +certain countries either by patents or by copyrighted interfaces, the +original copyright holder who places the Program under this License +may add an explicit geographical distribution limitation excluding +those countries, so that distribution is permitted only in or among +countries not thus excluded. In such case, this License incorporates +the limitation as if written in the body of this License. + + 9. The Free Software Foundation may publish revised and/or new versions +of the General Public License from time to time. Such new versions will +be similar in spirit to the present version, but may differ in detail to +address new problems or concerns. + +Each version is given a distinguishing version number. If the Program +specifies a version number of this License which applies to it and "any +later version", you have the option of following the terms and conditions +either of that version or of any later version published by the Free +Software Foundation. If the Program does not specify a version number of +this License, you may choose any version ever published by the Free Software +Foundation. + + 10. If you wish to incorporate parts of the Program into other free +programs whose distribution conditions are different, write to the author +to ask for permission. For software which is copyrighted by the Free +Software Foundation, write to the Free Software Foundation; we sometimes +make exceptions for this. Our decision will be guided by the two goals +of preserving the free status of all derivatives of our free software and +of promoting the sharing and reuse of software generally. + + NO WARRANTY + + 11. BECAUSE THE PROGRAM IS LICENSED FREE OF CHARGE, THERE IS NO WARRANTY +FOR THE PROGRAM, TO THE EXTENT PERMITTED BY APPLICABLE LAW. EXCEPT WHEN +OTHERWISE STATED IN WRITING THE COPYRIGHT HOLDERS AND/OR OTHER PARTIES +PROVIDE THE PROGRAM "AS IS" WITHOUT WARRANTY OF ANY KIND, EITHER EXPRESSED +OR IMPLIED, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF +MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE. THE ENTIRE RISK AS +TO THE QUALITY AND PERFORMANCE OF THE PROGRAM IS WITH YOU. SHOULD THE +PROGRAM PROVE DEFECTIVE, YOU ASSUME THE COST OF ALL NECESSARY SERVICING, +REPAIR OR CORRECTION. + + 12. IN NO EVENT UNLESS REQUIRED BY APPLICABLE LAW OR AGREED TO IN WRITING +WILL ANY COPYRIGHT HOLDER, OR ANY OTHER PARTY WHO MAY MODIFY AND/OR +REDISTRIBUTE THE PROGRAM AS PERMITTED ABOVE, BE LIABLE TO YOU FOR DAMAGES, +INCLUDING ANY GENERAL, SPECIAL, INCIDENTAL OR CONSEQUENTIAL DAMAGES ARISING +OUT OF THE USE OR INABILITY TO USE THE PROGRAM (INCLUDING BUT NOT LIMITED +TO LOSS OF DATA OR DATA BEING RENDERED INACCURATE OR LOSSES SUSTAINED BY +YOU OR THIRD PARTIES OR A FAILURE OF THE PROGRAM TO OPERATE WITH ANY OTHER +PROGRAMS), EVEN IF SUCH HOLDER OR OTHER PARTY HAS BEEN ADVISED OF THE +POSSIBILITY OF SUCH DAMAGES. + + END OF TERMS AND CONDITIONS + + How to Apply These Terms to Your New Programs + + If you develop a new program, and you want it to be of the greatest +possible use to the public, the best way to achieve this is to make it +free software which everyone can redistribute and change under these terms. + + To do so, attach the following notices to the program. It is safest +to attach them to the start of each source file to most effectively +convey the exclusion of warranty; and each file should have at least +the "copyright" line and a pointer to where the full notice is found. + + + Copyright (C) + + 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 2 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, write to the Free Software Foundation, Inc., + 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + +Also add information on how to contact you by electronic and paper mail. + +If the program is interactive, make it output a short notice like this +when it starts in an interactive mode: + + Gnomovision version 69, Copyright (C) year name of author + Gnomovision comes with ABSOLUTELY NO WARRANTY; for details type `show w'. + This is free software, and you are welcome to redistribute it + under certain conditions; type `show c' for details. + +The hypothetical commands `show w' and `show c' should show the appropriate +parts of the General Public License. Of course, the commands you use may +be called something other than `show w' and `show c'; they could even be +mouse-clicks or menu items--whatever suits your program. + +You should also get your employer (if you work as a programmer) or your +school, if any, to sign a "copyright disclaimer" for the program, if +necessary. Here is a sample; alter the names: + + Yoyodyne, Inc., hereby disclaims all copyright interest in the program + `Gnomovision' (which makes passes at compilers) written by James Hacker. + + , 1 April 1989 + Ty Coon, President of Vice + +This General Public License does not permit incorporating your program into +proprietary programs. If your program is a subroutine library, you may +consider it more useful to permit linking proprietary applications with the +library. If this is what you want to do, use the GNU Lesser General +Public License instead of this License. diff --git a/lib/__init__.py b/lib/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/lib/__init__.pyc b/lib/__init__.pyc new file mode 100644 index 0000000..9334b0b Binary files /dev/null and b/lib/__init__.pyc differ diff --git a/lib/analysetube.py b/lib/analysetube.py new file mode 100644 index 0000000..a60e63f --- /dev/null +++ b/lib/analysetube.py @@ -0,0 +1,175 @@ +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +# +## +# +#since v0.073 all cavities will start at (choosen axis) = 1 (-> cavity_coord) +class AnalyseTube(): + __groups = [] + __chain = () + __voxelPlanes = [] + def __init__ (self, groups, chain): + self.__chain = int(chain) #### chain ? ? ?? ?? this should be cavity !!! + self.__groups = groups + # self.__method ='' + def analyseTubeArea(self, axis): + print "Scanning cavity " + str(self.__chain)+" in "+axis+" direction areasize" + if (axis == 'x'): + self.__voxelPlanes = [0 for m in range (len (self.__groups))] + for xCount in range (len (self.__groups)): + for yCount in range (len (self.__groups[0])): + for zCount in range (len (self.__groups[0][0])): + if str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) != -1: + self.__voxelPlanes [xCount] += 1 + if (axis == 'y'): + self.__voxelPlanes = [0 for m in range (len (self.__groups[0]))] + for xCount in range (len (self.__groups)): + for yCount in range (len (self.__groups[0])): + for zCount in range (len (self.__groups[0][0])): + if str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) != -1: + self.__voxelPlanes [yCount] += 1 + if (axis == 'z'): + self.__voxelPlanes = [0 for m in range (len (self.__groups[0][0]))] + for xCount in range (len (self.__groups)): + for yCount in range (len (self.__groups[0])): + for zCount in range (len (self.__groups[0][0])): + if str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) != -1: + self.__voxelPlanes [zCount] += 1 + + print self.__voxelPlanes + def analyseTubeDiameterMax(self, axis): + print "Scanning cavity " + str(self.__chain)+" in "+axis+" direction for maximum diameter" + + if (axis == 'x'): + self.__voxelPlanes = [0 for m in range (len (self.__groups))] + for xCount in range (len (self.__groups)): + maximum = 0 +#find maximum diameter in z direction + for yCount in range (len (self.__groups[0])): + sum = 0 + for zCount in range (len (self.__groups[0][0])): +#if a gap in the tunnel is found espacily at the corners ... + if (sum > 0 and str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) == -1): + if sum > maximum: + maximum = sum + sum = 0 + if str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) != -1: + sum += 1 + if sum > maximum: + maximum = sum + self.__voxelPlanes [xCount] = maximum +#find maximum diameter in y direction + for zCount in range (len (self.__groups[0][0])): + sum = 0 + for yCount in range (len (self.__groups[0])): +#if a gap in the tunnel is found espacily at the corners ... + if (sum > 0 and str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) == -1): + if sum > maximum: + maximum = sum + sum = 0 + if str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) != -1: + sum += 1 + if sum > maximum: + maximum = sum + self.__voxelPlanes [xCount] = maximum + + if (axis == 'y'): + self.__voxelPlanes = [0 for m in range (len (self.__groups[0]))] + for yCount in range (len (self.__groups[0])): + maximum = 0 +#find maximum diameter in z direction + for xCount in range (len (self.__groups)): + sum = 0 + for zCount in range (len (self.__groups[0][0])): +#if a gap in the tunnel is found espacily at the corners ... + if (sum > 0 and str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) == -1): + if sum > maximum: + maximum = sum + sum = 0 + if str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) != -1: + sum += 1 + if sum > maximum: + maximum = sum + self.__voxelPlanes [yCount] = maximum +#find maximum diameter in x direction + for zCount in range (len (self.__groups[0][0])): + sum = 0 + for xCount in range (len (self.__groups)): +#if a gap in the tunnel is found espacily at the corners ... + if (sum > 0 and str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) == -1): + if sum > maximum: + maximum = sum + sum = 0 + if str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) != -1: + sum += 1 + if sum > maximum: + maximum = sum + self.__voxelPlanes [yCount] = maximum + + if (axis == 'z'): + self.__voxelPlanes = [0 for m in range (len (self.__groups[0][0]))] + for zCount in range (len (self.__groups[0][0])): + maximum = 0 +#find maximum diameter in x direction + for yCount in range (len (self.__groups[0])): + sum = 0 + for xCount in range (len (self.__groups)): +#if a gap in the tunnel is found espacily at the corners ... + #if (sum > 0 and self.__groups[xCount-1][yCount][zCount] != self.__chain): + if (sum > 0 and str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) == -1): + if sum > maximum: + maximum = sum + sum = 0 + if str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) != -1: + sum += 1 + if sum > maximum: + maximum = sum + self.__voxelPlanes [zCount] = maximum + +#find maximum diameter in y direction + for xCount in range (len (self.__groups)): + sum = 0 + for yCount in range (len (self.__groups[0])): +#if a gap in the tunnel is found espacily at the corners ... + if (sum > 0 and str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) == -1): + if sum > maximum: + maximum = sum + sum = 0 + if str(self.__groups[xCount][yCount][zCount]).find(str(self.__chain)) != -1: + sum += 1 + if sum > maximum: + maximum = sum + self.__voxelPlanes [zCount] = maximum + print self.__voxelPlanes + +# since v0.073 all cavities will start at (choosen axis) = 1 (-> cavity_coord) + def writeStatistics(self, file): + statisticFH = open (file,'w') + cavity_coord = 0 + for axisCount in range (len (self.__voxelPlanes)): + if cavity_coord > 0: + cavity_coord += 1 + if self.__voxelPlanes[axisCount] > 0: + if cavity_coord == 0: + cavity_coord = 1 + statisticFH.write (str(cavity_coord)+" " +str(self.__voxelPlanes[axisCount])+"\n") +#since + else: + statisticFH.write (str(cavity_coord)+" 0\n") + statisticFH.close() + print "Saved analysed cavity "+ str(self.__chain) +" in\n" + file + diff --git a/lib/analysetube.pyc b/lib/analysetube.pyc new file mode 100644 index 0000000..7faa802 Binary files /dev/null and b/lib/analysetube.pyc differ diff --git a/lib/averagedx.py b/lib/averagedx.py new file mode 100644 index 0000000..5aa47a2 --- /dev/null +++ b/lib/averagedx.py @@ -0,0 +1,102 @@ +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +# +## + + + +#### +# not implemented yet +# first draft +# untestet +### +from buildopendx import * +class averageDx(): + __origin = [] + __densityObjectArray = [] + __dimention = [] + __average = [] + def __init__(self, dxObjectArray): + self.__densityObjectArray = dxObjectArray + self.calcNewOrigin() + self.calcNewDimention() + self.__calcAverages() +# getting the lowest origin : + def calcNewOrigin(self): + origin = self.__densityObjectArray[0].getOrigin() + for i in range(len (self.__densityObjectArray)): +#x + if (self.__densityObjectArray[i].getOrigin()[0] < origin[0]): + origin[0] = self.__densityObjectArray[i].getOrigin()[0] +#y + if (self.__densityObjectArray[i].getOrigin()[1] < origin[1]): + origin[1] = self.__densityObjectArray[i].getOrigin()[1] +#z + if (self.__densityObjectArray[i].getOrigin()[2] < origin[0]): + origin[2] = self.__densityObjectArray[i].getOrigin()[2] + self.__origin = origin + print "New Origin:\n"+str(origin[0])+" "+str(origin[1])+" "+str(origin[2]) + +#calc new size: + def calcNewDimention(self): + dimention = self.__densityObjectArray[0].getDimention() + for i in range(len (self.__densityObjectArray)): +#x + if (self.__densityObjectArray[i].getDimention()[0] < dimention[0]): + dimention[0] = self.__densityObjectArray[i].getDimention()[0] +#y + if (self.__densityObjectArray[i].getDimention()[1] < Dimention[1]): + dimention[1] = self.__densityObjectArray[i].getDimention()[1] +#z + if (self.__densityObjectArray[i].getDimention()[2] < Dimention[0]): + dimention[2] = self.__densityObjectArray[i].getDimention()[2] + self.__dimention = dimention + print "New Dimention:\n"+str(dimention[0])+" "+str(dimention[1])+" "+str(dimention[2]) + def __calcAverages(self): + sum =[[[0 for zR in range(self.__origin[0][0])] for yR in range(self.__origin[0])] for xR in range(self.__origin)] + average =[[[0 for zR in range(self.__origin[0][0])] for yR in range(self.__origin[0])] for xR in range(self.__origin)] +#looping over all dx files + for dxObjectCount in range(len (self.__densityObjectArray)): +#loop over the new grid + print "Start summing densities on file"+ str(dxObjectCount) + for xCount in range (len (self.__dimention)): + for yCount in range (len (self.__dimention[0])): + for zCount in range (len (self.__dimention[0][0])): +#testing if current voxel lies inside of current dx grid +# "left side" + if (self.__origin[0] + xCount > self.__densityObjectArray[dxObjectCount].getOrigin()[0] and + self.__origin[1] + yCount > self.__densityObjectArray[dxObjectCount].getOrigin()[1] and + self.__origin[2] + zCount > self.__densityObjectArray[dxObjectCount].getOrigin()[2] and +#"right side" + self.__origin[0] + xCount < self.__densityObjectArray[dxObjectCount].getOrigin()[0] + self.__densityObjectArray[dxObjectCount].getDimention[0] and + self.__origin[1] + yCount < self.__densityObjectArray[dxObjectCount].getOrigin()[1] + self.__densityObjectArray[dxObjectCount].getDimention[1] and + self.__origin[2] + zCount < self.__densityObjectArray[dxObjectCount].getOrigin()[2] + self.__densityObjectArray[dxObjectCount].getDimention[2] ): +#calculating shift in x y z direction + xShift = self.__origin[0] - self.__densityObjectArray[dxObjectCount].getOrigin()[0] + yShift = self.__origin[1] - self.__densityObjectArray[dxObjectCount].getOrigin()[1] + zShift = self.__origin[2] - self.__densityObjectArray[dxObjectCount].getOrigin()[2] + sum[xCount][yCount][zCount] += self.__densityObjectArray[dxObjectCount].getDensity()[xShift+xCount][yShift+yCount][zShift+zCount] + print "Calculating averages" + for xCount in range (len (self.__dimention)): + for yCount in range (len (self.__dimention[0])): + for zCount in range (len (self.__dimention[0][0])): + average [xCount][yCount][zCount]=sum[xCount][yCount][zCount] / (len (self.__densityObjectArray) + self.__average = average + def getDxObject(self=None): + object = BuildOpenDx (self.__average, self.__densityObjectArray[0].getStepsize(), self.__origin, self.__dimention) + return object + diff --git a/lib/buildopendx.py b/lib/buildopendx.py new file mode 100644 index 0000000..79ab3a8 --- /dev/null +++ b/lib/buildopendx.py @@ -0,0 +1,50 @@ +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +# +## + +##### +# v0.15 +# +setDensity +# +#since v0.065 dimention is a new mandatory parameter! +##### +class BuildOpenDx(): +#density array [x][y][z] = density as float [dx format] + __density=[] +#origin array [0] = x, [1] =y, [2] = z as float + __origin =[] + __stepSize = () +#dimention array [0] = x, [1] =y, [2] = z as integer + __dimention = [] + def __init__ (self, density, stepSize, origin,dimention): + self.__density = density + self.__stepSize = stepSize + self.__origin = origin + self.__dimention = dimention + + def getOrigin(self = None): + return self.__origin + def getDensity(self=None): + return self.__density + def setDensity(self, density): + self.__density = density + + def getStepsize(self=None): + return self.__stepSize + def getDimention(self=None): + return self.__dimention diff --git a/lib/buildopendx.pyc b/lib/buildopendx.pyc new file mode 100644 index 0000000..7cf3b33 Binary files /dev/null and b/lib/buildopendx.pyc differ diff --git a/lib/filltube.py b/lib/filltube.py new file mode 100644 index 0000000..4360c2b --- /dev/null +++ b/lib/filltube.py @@ -0,0 +1,724 @@ +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + + +from buildopendx import * +####### +# This class fills cavities along principle axis after grouping +# +# v0.12 +# +int(scanMethod) +######## +class FillTube(): + __dxObjectTube =() + __dxObjectProtein =() + __scanMethod =() + __minProteinDensity =() + __newTubeDensity = [] +#voxelGroupGrid +#[x][y][z] = group + __voxelGroupGrid =() +#groups +#[[x][y][z][x][y][z]] + __groups =() + def __init__(self,dxObjectWater,dxObjectProtein, dxObjectTube, scanMethod, minProteinDensity, minSolvDensity, groups = None, voxelGroupGrid = None): + self.__dxObjectTube = dxObjectTube + self.__dxObjectWater = dxObjectWater + self.__dxObjectProtein = dxObjectProtein + self.__scanMethod = int (scanMethod) + self.__minProteinDensity = float(minProteinDensity) + self.__groups = groups + self.__voxelGroupGrid = voxelGroupGrid + self.__minSolvDensity = float(minSolvDensity) + def getNewTubeObject(self): + object = BuildOpenDx (self.__newTubeDensity, self.__dxObjectTube.getStepsize(), self.__dxObjectTube.getOrigin(), self.__dxObjectTube.getDimention()) + return object + def getVoxelGroupGrid(self): + return self.__voxelGroupGrid + def getGroups(self): + return self.__groups + + + def fillGroups(self): + if self.__scanMethod == 2: + self.__fill2d() + if self.__scanMethod == 3: + self.__fill3d() + def __fill2d(self): + print "Filling cavities after 2D scanning" + oldTube = self.__dxObjectTube.getDensity() + newTube = [[[-1 for zR in range(len(oldTube[0][0]))] for yR in range(len(oldTube[0]))] for xR in range(len(oldTube))] + newVoxelGroupGrid = [[[-1 for zR in range(len(oldTube[0][0]))] for yR in range(len(oldTube[0]))] for xR in range(len(oldTube))] + densityProtein = self.__dxObjectProtein.getDensity() + densityWater = self.__dxObjectWater.getDensity() +############ X ############### + print "X" + for yCount in range (len (densityProtein[0])): + for zCount in range (len (densityProtein[0][0])): + for xCount in range (len (densityProtein)): +#### filling X+ +#### xPlus xMinus indicates if a cavity will be filled in plus or minus direction + + xPlus = 1 + xMinus = 1 + if xCount == (len (densityProtein)-1): + xPlus = 0 + elif (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity) or + densityProtein[xCount+1][yCount][zCount] > float(self.__minProteinDensity) ): + xPlus = 0 + +####### this block is taken from planesearch to indicate towards which principle directions cavities can be filled + xProteinLeft = 0 #x- + xProteinRight = 0 #x+ + j = xCount + while (j > 0): + j -=1 + if (densityProtein[j][yCount][zCount] >= float(self.__minProteinDensity)): + xProteinLeft = xCount - j + break + + for k in range (xCount,len(densityProtein)): + if (densityProtein[k][yCount][zCount] >= float(self.__minProteinDensity) ): + xProteinRight = k- xCount + break + if xProteinLeft < 1: + xMinus = 0 + if xProteinRight <1: + xPlus = 0 +#### end of planesearch block#### + + + +#calculating shift between water & tube coordienates towards protein coordinates + proteinOrigin = self.__dxObjectProtein.getOrigin() + tubeOrigin = self.__dxObjectTube.getOrigin() + tubeToProteinOrigin = (int (proteinOrigin[0] - tubeOrigin[0]),int (proteinOrigin[1] - tubeOrigin[1]),int (proteinOrigin[2] - tubeOrigin[2])) + convertPosition = (xCount+tubeToProteinOrigin[0], yCount+tubeToProteinOrigin[1],zCount+tubeToProteinOrigin[2]) + if (oldTube[convertPosition[0]+1][convertPosition[1]][convertPosition[2]] > 0): + xPlus = 0 + if (densityWater[convertPosition[0]+1][convertPosition[1]][convertPosition[2]] < self.__minSolvDensity): + xPlus = 0 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]] > 0 and xPlus == 1): + i = 1 + oldGroup = self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]] + while ( i > -1): + if (densityWater[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] > self.__minSolvDensity): + newTube[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = densityWater[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] +#adding new voxel to existing group definitions + if self.__voxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] < 0 or self.__voxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] == oldGroup: + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = oldGroup + else: + tmp = str(self.__voxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]]) + tmp += " "+str(oldGroup) + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = tmp +## adding voxel to groups +#[[x][y][z][x][y][z]] + x = convertPosition[0]+i + y = convertPosition[1] + z = convertPosition[2] + self.__groups[oldGroup].append([x,y,z]) + i +=1 + if (oldTube[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] > 0): + i = -1 + if (densityWater[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] < self.__minSolvDensity): + i = -1 + if (densityProtein[xCount+i][yCount][zCount] > float(self.__minProteinDensity)): + i = -1 + +#### filling X- + + if xCount == 0: + xMinus = 0 + elif (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity) or + #if (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity) or + densityProtein[xCount-1][yCount][zCount] > float(self.__minProteinDensity) ): + xMinus = 0 +#calculating shift between water & tube coordienates towards protein coordinates + proteinOrigin = self.__dxObjectProtein.getOrigin() + tubeOrigin = self.__dxObjectTube.getOrigin() + tubeToProteinOrigin = (int (proteinOrigin[0] - tubeOrigin[0]),int (proteinOrigin[1] - tubeOrigin[1]),int (proteinOrigin[2] - tubeOrigin[2])) + convertPosition = (xCount+tubeToProteinOrigin[0], yCount+tubeToProteinOrigin[1],zCount+tubeToProteinOrigin[2]) + if (oldTube[convertPosition[0]-1][convertPosition[1]][convertPosition[2]] > 0): + xMinus = 0 + if (densityWater[convertPosition[0]-1][convertPosition[1]][convertPosition[2]] < self.__minSolvDensity): + xMinus = 0 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]] > 0 and xMinus == 1): + i = -1 + oldGroup = self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]] + while ( i < 1): + if (densityWater[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] > self.__minSolvDensity): + newTube[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = densityWater[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] +#adding new voxel to existing group definitions + if self.__voxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] < 0 or self.__voxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] == oldGroup: + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = oldGroup + else: + tmp = str(self.__voxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]]) + tmp += " "+str(oldGroup) + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = tmp +## adding voxel to groups +#[[x][y][z][x][y][z]] + x = convertPosition[0]+i + y = convertPosition[1] + z = convertPosition[2] + self.__groups[oldGroup].append([x,y,z]) + i -=1 + if (oldTube[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] > 0): + i = 1 + if (densityWater[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] < self.__minSolvDensity): + i = 1 + if (densityProtein[xCount+i][yCount][zCount] > float(self.__minProteinDensity)): + i = 1 + + + + + +############ Y ############### + print "Y" + for xCount in range (len (densityProtein)): + for zCount in range (len (densityProtein[0][0])): + for yCount in range (len (densityProtein[0])): + +##### filling Y+ + yPlus = 1 + yMinus = 1 + + yProteinLeft = 0 #y- + yProteinRight = 0 #y+ + j = yCount + while (j > 0): + j -=1 + if (densityProtein[xCount][j][zCount] >= float(self.__minProteinDensity) ): + yProteinLeft = yCount - j + break + + for k in range (yCount,len(densityProtein[0])): + if (densityProtein[xCount][k][zCount] >= float(self.__minProteinDensity) ): + yProteinRight = k - yCount + break + if yProteinLeft < 1: + yMinus = 0 + if xProteinRight <1: + yPlus = 0 + + + + if yCount == (len (densityProtein[0])-1): + yPlus = 0 +# if (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity)or +# densityProtein[xCount][yCount+1][zCount] > float(self.__minProteinDensity) ): + elif (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity)or + densityProtein[xCount][yCount+1][zCount] > float(self.__minProteinDensity) ): + yPlus = 0 + +#calculating shift between water & tube coordienates towards protein coordinates + proteinOrigin = self.__dxObjectProtein.getOrigin() + tubeOrigin = self.__dxObjectTube.getOrigin() + tubeToProteinOrigin = (int (proteinOrigin[0] - tubeOrigin[0]),int (proteinOrigin[1] - tubeOrigin[1]),int (proteinOrigin[2] - tubeOrigin[2])) + convertPosition = (xCount+tubeToProteinOrigin[0], yCount+tubeToProteinOrigin[1],zCount+tubeToProteinOrigin[2]) + if (oldTube[convertPosition[0]][convertPosition[1]+1][convertPosition[2]] > 0): + yPlus = 0 + if (densityWater[convertPosition[0]][convertPosition[1]+1][convertPosition[2]] < self.__minSolvDensity): + yPlus = 0 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]] > 0 and yPlus ==1): + i = 1 + oldGroup = self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]] + while ( i > -1): + if (densityWater[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] > self.__minSolvDensity): + newTube[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] = densityWater[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] +#adding new voxel to existing group definitions + if self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] < 0 or self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] == oldGroup: + newVoxelGroupGrid[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] = oldGroup + else: + tmp = str(self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i]) + tmp += " "+str(oldGroup) + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = tmp +## adding voxel to groups +#[[x][y][z][x][y][z]] + x = convertPosition[0] + y = convertPosition[1]+i + z = convertPosition[2] + self.__groups[oldGroup].append([x,y,z]) + i +=1 + if (densityWater[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] < self.__minSolvDensity): + i = -1 + if (oldTube[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] > 0): + i = -1 + if (densityProtein[xCount][yCount+i][zCount] > float(self.__minProteinDensity)): + i = -1 + + +###### filling Y- + + if yCount == 0: + yMinus = 0 + elif (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity)or + densityProtein[xCount][yCount-1][zCount] > float(self.__minProteinDensity) ): + yMinus = 0 + +#calculating shift between water & tube coordienates towards protein coordinates + proteinOrigin = self.__dxObjectProtein.getOrigin() + tubeOrigin = self.__dxObjectTube.getOrigin() + tubeToProteinOrigin = (int (proteinOrigin[0] - tubeOrigin[0]),int (proteinOrigin[1] - tubeOrigin[1]),int (proteinOrigin[2] - tubeOrigin[2])) + convertPosition = (xCount+tubeToProteinOrigin[0], yCount+tubeToProteinOrigin[1],zCount+tubeToProteinOrigin[2]) + if (oldTube[convertPosition[0]][convertPosition[1]-1][convertPosition[2]] > 0): + yMinus = 0 + if (densityWater[convertPosition[0]][convertPosition[1]-1][convertPosition[2]] < self.__minSolvDensity): + yMinus = 0 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]] > 0 and yMinus == 1): + i = -1 + oldGroup = self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]] + while ( i < 1): + if (densityWater[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] > self.__minSolvDensity): + newTube[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] = densityWater[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] +#adding new voxel to existing group definitions + if self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] < 0 or self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] == oldGroup: + newVoxelGroupGrid[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] = oldGroup + else: + tmp = str(self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i]) + tmp += " "+str(oldGroup) + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = tmp +## adding voxel to groups +#[[x][y][z][x][y][z]] + x = convertPosition[0] + y = convertPosition[1]+i + z = convertPosition[2] + self.__groups[oldGroup].append([x,y,z]) + i -=1 + if (densityWater[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] < self.__minSolvDensity): + i = 1 + if (oldTube[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] > 0): + i = 1 + if (densityProtein[xCount][yCount+i][zCount] > float(self.__minProteinDensity)): + i = 1 + + +####### Z ######## + print "Z" + for xCount in range (len (densityProtein)): + for yCount in range (len (densityProtein[0])): + for zCount in range (len (densityProtein[0][0])): + +########## filling Z+ + zPlus = 1 + zMinus = 1 + + zProteinLeft = 0#z- + zProteinRight = 0#z+ + j = zCount + while (j > 0): + j -=1 + if (densityProtein[xCount][yCount][j] >= float(self.__minProteinDensity) ): + zProteinLeft = zCount - j + break + + for k in range (zCount,len(densityProtein[0][0])): + if (densityProtein[xCount][yCount][k] >= float(self.__minProteinDensity) ): + zProteinRight =k- zCount + break + + if zProteinLeft < 1: + zMinus = 0 + if zProteinRight <1: + zPlus = 0 + + if zCount == (len (densityProtein[0][0])-1): + zPlus = 0 + elif (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity) or + densityProtein[xCount][yCount][zCount+1] > float(self.__minProteinDensity) ): + zPlus = 0 +#calculating shift between water & tube coordienates towards protein coordinates + proteinOrigin = self.__dxObjectProtein.getOrigin() + tubeOrigin = self.__dxObjectTube.getOrigin() + tubeToProteinOrigin = (int (proteinOrigin[0] - tubeOrigin[0]),int (proteinOrigin[1] - tubeOrigin[1]),int (proteinOrigin[2] - tubeOrigin[2])) + convertPosition = (xCount+tubeToProteinOrigin[0], yCount+tubeToProteinOrigin[1],zCount+tubeToProteinOrigin[2]) + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]+1] > 0): + zPlus = 0 + if (densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]+1] < self.__minSolvDensity): + zPlus = 0 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]] > 0 and zPlus == 1): + i = 1 + oldGroup = self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]] + while ( i > -1): + if (densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] > self.__minSolvDensity): + newTube[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] = densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] +#adding new voxel to existing group definitions + if self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] < 0 or self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] == oldGroup: + newVoxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] = oldGroup + else: + tmp = str(self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i]) + tmp += " "+str(oldGroup) + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = tmp +## adding voxel to groups +#[[x][y][z][x][y][z]] + x = convertPosition[0] + y = convertPosition[1] + z = convertPosition[2]+i + self.__groups[oldGroup].append([x,y,z]) + i +=1 + if (densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] < self.__minSolvDensity): + i = -1 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] > 0): + i = -1 + if (densityProtein[xCount][yCount][zCount+i] > float(self.__minProteinDensity)): + i = -1 +########### filliung Z- + + if zCount == 0: + zMinus = 0 + elif (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity) or + densityProtein[xCount][yCount][zCount-1] > float(self.__minProteinDensity) ): + zMinus = 0 +#calculating shift between water & tube coordienates towards protein coordinates + proteinOrigin = self.__dxObjectProtein.getOrigin() + tubeOrigin = self.__dxObjectTube.getOrigin() + tubeToProteinOrigin = (int (proteinOrigin[0] - tubeOrigin[0]),int (proteinOrigin[1] - tubeOrigin[1]),int (proteinOrigin[2] - tubeOrigin[2])) + convertPosition = (xCount+tubeToProteinOrigin[0], yCount+tubeToProteinOrigin[1],zCount+tubeToProteinOrigin[2]) + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]-1] > 0): + zMinus = 0 + if (densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]-1] < self.__minSolvDensity): + zMinus = 0 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]] > 0 and zMinus == 1): + i = -1 + while ( i < 0): + if (densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] > self.__minSolvDensity): + newTube[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] = densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] + if self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] < 0 or self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] == oldGroup: + newVoxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] = oldGroup + else: + tmp = str(self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i]) + tmp += " "+str(oldGroup) + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = tmp + ## adding voxel to groups + #[[x][y][z][x][y][z]] + x = convertPosition[0] + y = convertPosition[1] + z = convertPosition[2]+i + self.__groups[oldGroup].append([x,y,z]) + i -=1 + if (densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] < self.__minSolvDensity): + i = 1 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] > 0): + i = 1 + if (densityProtein[xCount][yCount][zCount+i] > float(self.__minProteinDensity)): + i = 1 + + + for xCount in range (len (newTube)): + for yCount in range (len (newTube[0])): + for zCount in range (len (newTube[0][0])): + if newTube[xCount][yCount][zCount] > 0: + oldTube[xCount][yCount][zCount] = newTube[xCount][yCount][zCount] + self.__voxelGroupGrid[xCount][yCount][zCount] = newVoxelGroupGrid[xCount][yCount][zCount] + self.__newTubeDensity = oldTube + + + def __fill3d(self): + print "Filling cavities after 3D scanning" + oldTube = self.__dxObjectTube.getDensity() + newTube = [[[-1 for zR in range(len(oldTube[0][0]))] for yR in range(len(oldTube[0]))] for xR in range(len(oldTube))] + newVoxelGroupGrid = [[[-1 for zR in range(len(oldTube[0][0]))] for yR in range(len(oldTube[0]))] for xR in range(len(oldTube))] + densityProtein = self.__dxObjectProtein.getDensity() + densityWater = self.__dxObjectWater.getDensity() +############ X ############### + print "X" + for yCount in range (len (densityProtein[0])): + for zCount in range (len (densityProtein[0][0])): + for xCount in range (len (densityProtein)): +#### filling X+ +#### xPlus xMinus indicates if a cavity will be filled in plus or minus direction + + xPlus = 1 + if xCount == (len (densityProtein)-1): + xPlus = 0 + elif (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity) or + densityProtein[xCount+1][yCount][zCount] > float(self.__minProteinDensity) ): + xPlus = 0 + +#calculating shift between water & tube coordienates towards protein coordinates + proteinOrigin = self.__dxObjectProtein.getOrigin() + tubeOrigin = self.__dxObjectTube.getOrigin() + tubeToProteinOrigin = (int (proteinOrigin[0] - tubeOrigin[0]),int (proteinOrigin[1] - tubeOrigin[1]),int (proteinOrigin[2] - tubeOrigin[2])) + convertPosition = (xCount+tubeToProteinOrigin[0], yCount+tubeToProteinOrigin[1],zCount+tubeToProteinOrigin[2]) + if (oldTube[convertPosition[0]+1][convertPosition[1]][convertPosition[2]] > 0): + xPlus = 0 + if (densityWater[convertPosition[0]+1][convertPosition[1]][convertPosition[2]] < self.__minSolvDensity): + xPlus = 0 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]] > 0 and xPlus == 1): + i = 1 + oldGroup = self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]] + while ( i > -1): + if (densityWater[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] > self.__minSolvDensity): + newTube[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = densityWater[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] +#adding new voxel to existing group definitions + if self.__voxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] < 0 or self.__voxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] == oldGroup: + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = oldGroup + else: + tmp = str(self.__voxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]]) + tmp += " "+str(oldGroup) + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = tmp +## adding voxel to groups +#[[x][y][z][x][y][z]] + x = convertPosition[0]+i + y = convertPosition[1] + z = convertPosition[2] + self.__groups[oldGroup].append([x,y,z]) + i +=1 + if (oldTube[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] > 0): + i = -1 + if (densityWater[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] < self.__minSolvDensity): + i = -1 + if (densityProtein[xCount+i][yCount][zCount] > float(self.__minProteinDensity)): + i = -1 + +#### filling X- + xMinus = 1 + if xCount == 0: + xMinus = 0 + elif (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity) or + #if (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity) or + densityProtein[xCount-1][yCount][zCount] > float(self.__minProteinDensity) ): + xMinus = 0 +#calculating shift between water & tube coordienates towards protein coordinates + proteinOrigin = self.__dxObjectProtein.getOrigin() + tubeOrigin = self.__dxObjectTube.getOrigin() + tubeToProteinOrigin = (int (proteinOrigin[0] - tubeOrigin[0]),int (proteinOrigin[1] - tubeOrigin[1]),int (proteinOrigin[2] - tubeOrigin[2])) + convertPosition = (xCount+tubeToProteinOrigin[0], yCount+tubeToProteinOrigin[1],zCount+tubeToProteinOrigin[2]) + if (oldTube[convertPosition[0]-1][convertPosition[1]][convertPosition[2]] > 0): + xMinus = 0 + if (densityWater[convertPosition[0]-1][convertPosition[1]][convertPosition[2]] < self.__minSolvDensity): + xMinus = 0 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]] > 0 and xMinus == 1): + i = -1 + oldGroup = self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]] + while ( i < 1): + if (densityWater[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] > self.__minSolvDensity): + newTube[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = densityWater[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] +#adding new voxel to existing group definitions + if self.__voxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] < 0 or self.__voxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] == oldGroup: + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = oldGroup + else: + tmp = str(self.__voxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]]) + tmp += " "+str(oldGroup) + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = tmp +## adding voxel to groups +#[[x][y][z][x][y][z]] + x = convertPosition[0]+i + y = convertPosition[1] + z = convertPosition[2] + self.__groups[oldGroup].append([x,y,z]) + i -=1 + if (oldTube[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] > 0): + i = 1 + if (densityWater[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] < self.__minSolvDensity): + i = 1 + if (densityProtein[xCount+i][yCount][zCount] > float(self.__minProteinDensity)): + i = 1 + + +############ Y ############### + print "Y" + for xCount in range (len (densityProtein)): + for zCount in range (len (densityProtein[0][0])): + for yCount in range (len (densityProtein[0])): + +##### filling Y+ + yPlus = 1 + if yCount == (len (densityProtein[0])-1): + yPlus = 0 +# if (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity)or +# densityProtein[xCount][yCount+1][zCount] > float(self.__minProteinDensity) ): + elif (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity)or + densityProtein[xCount][yCount+1][zCount] > float(self.__minProteinDensity) ): + yPlus = 0 + +#calculating shift between water & tube coordienates towards protein coordinates + proteinOrigin = self.__dxObjectProtein.getOrigin() + tubeOrigin = self.__dxObjectTube.getOrigin() + tubeToProteinOrigin = (int (proteinOrigin[0] - tubeOrigin[0]),int (proteinOrigin[1] - tubeOrigin[1]),int (proteinOrigin[2] - tubeOrigin[2])) + convertPosition = (xCount+tubeToProteinOrigin[0], yCount+tubeToProteinOrigin[1],zCount+tubeToProteinOrigin[2]) + if (oldTube[convertPosition[0]][convertPosition[1]+1][convertPosition[2]] > 0): + yPlus = 0 + if (densityWater[convertPosition[0]][convertPosition[1]+1][convertPosition[2]] < self.__minSolvDensity): + yPlus = 0 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]] > 0 and yPlus ==1): + i = 1 + oldGroup = self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]] + while ( i > -1): + if (densityWater[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] > self.__minSolvDensity): + newTube[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] = densityWater[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] +#adding new voxel to existing group definitions + if self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] < 0 or self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] == oldGroup: + newVoxelGroupGrid[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] = oldGroup + else: + tmp = str(self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i]) + tmp += " "+str(oldGroup) + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = tmp +## adding voxel to groups +#[[x][y][z][x][y][z]] + x = convertPosition[0] + y = convertPosition[1]+i + z = convertPosition[2] + self.__groups[oldGroup].append([x,y,z]) + i +=1 + if (densityWater[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] < self.__minSolvDensity): + i = -1 + if (oldTube[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] > 0): + i = -1 + if (densityProtein[xCount][yCount+i][zCount] > float(self.__minProteinDensity)): + i = -1 + +###### filling Y- + yMinus = 1 + if yCount == 0: + yMinus = 0 + elif (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity)or + densityProtein[xCount][yCount-1][zCount] > float(self.__minProteinDensity) ): + yMinus = 0 + +#calculating shift between water & tube coordienates towards protein coordinates + proteinOrigin = self.__dxObjectProtein.getOrigin() + tubeOrigin = self.__dxObjectTube.getOrigin() + tubeToProteinOrigin = (int (proteinOrigin[0] - tubeOrigin[0]),int (proteinOrigin[1] - tubeOrigin[1]),int (proteinOrigin[2] - tubeOrigin[2])) + convertPosition = (xCount+tubeToProteinOrigin[0], yCount+tubeToProteinOrigin[1],zCount+tubeToProteinOrigin[2]) + if (oldTube[convertPosition[0]][convertPosition[1]-1][convertPosition[2]] > 0): + yMinus = 0 + if (densityWater[convertPosition[0]][convertPosition[1]-1][convertPosition[2]] < self.__minSolvDensity): + yMinus = 0 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]] > 0 and yMinus == 1): + i = -1 + oldGroup = self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]] + while ( i < 1): + if (densityWater[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] > self.__minSolvDensity): + newTube[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] = densityWater[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] +#adding new voxel to existing group definitions + if self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] < 0 or self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] == oldGroup: + newVoxelGroupGrid[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] = oldGroup + else: + tmp = str(self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i]) + tmp += " "+str(oldGroup) + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = tmp +## adding voxel to groups +#[[x][y][z][x][y][z]] + x = convertPosition[0] + y = convertPosition[1]+i + z = convertPosition[2] + self.__groups[oldGroup].append([x,y,z]) + i -=1 + if (densityWater[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] < self.__minSolvDensity): + i = 1 + if (oldTube[convertPosition[0]][convertPosition[1]+i][convertPosition[2]] > 0): + i = 1 + if (densityProtein[xCount][yCount+i][zCount] > float(self.__minProteinDensity)): + i = 1 + + +####### Z ######## + print "Z" + for xCount in range (len (densityProtein)): + for yCount in range (len (densityProtein[0])): + for zCount in range (len (densityProtein[0][0])): + +########## filling Z+ + zPlus = 1 + if zCount == (len (densityProtein[0][0])-1): + zPlus = 0 + elif (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity) or + densityProtein[xCount][yCount][zCount+1] > float(self.__minProteinDensity) ): + zPlus = 0 +#calculating shift between water & tube coordienates towards protein coordinates + proteinOrigin = self.__dxObjectProtein.getOrigin() + tubeOrigin = self.__dxObjectTube.getOrigin() + tubeToProteinOrigin = (int (proteinOrigin[0] - tubeOrigin[0]),int (proteinOrigin[1] - tubeOrigin[1]),int (proteinOrigin[2] - tubeOrigin[2])) + convertPosition = (xCount+tubeToProteinOrigin[0], yCount+tubeToProteinOrigin[1],zCount+tubeToProteinOrigin[2]) + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]+1] > 0): + zPlus = 0 + if (densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]+1] < self.__minSolvDensity): + zPlus = 0 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]] > 0 and zPlus == 1): + i = 1 + oldGroup = self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]] + while ( i > -1): + if (densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] > self.__minSolvDensity): + newTube[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] = densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] +#adding new voxel to existing group definitions + if self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] < 0 or self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] == oldGroup: + newVoxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] = oldGroup + else: + tmp = str(self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i]) + tmp += " "+str(oldGroup) + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = tmp +## adding voxel to groups +#[[x][y][z][x][y][z]] + x = convertPosition[0] + y = convertPosition[1] + z = convertPosition[2]+i + self.__groups[oldGroup].append([x,y,z]) + i +=1 + if (densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] < self.__minSolvDensity): + i = -1 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] > 0): + i = -1 + if (densityProtein[xCount][yCount][zCount+i] > float(self.__minProteinDensity)): + i = -1 +########### filliung Z- + zMinus = 1 + if zCount == 0: + zMinus = 0 + elif (densityProtein[xCount][yCount][zCount] > float(self.__minProteinDensity) or + densityProtein[xCount][yCount][zCount-1] > float(self.__minProteinDensity) ): + zMinus = 0 +#calculating shift between water & tube coordienates towards protein coordinates + proteinOrigin = self.__dxObjectProtein.getOrigin() + tubeOrigin = self.__dxObjectTube.getOrigin() + tubeToProteinOrigin = (int (proteinOrigin[0] - tubeOrigin[0]),int (proteinOrigin[1] - tubeOrigin[1]),int (proteinOrigin[2] - tubeOrigin[2])) + convertPosition = (xCount+tubeToProteinOrigin[0], yCount+tubeToProteinOrigin[1],zCount+tubeToProteinOrigin[2]) + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]-1] > 0): + zMinus = 0 + if (densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]-1] < self.__minSolvDensity): + zMinus = 0 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]] > 0 and zMinus == 1): + i = -1 + while ( i < 0): + if (densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] > self.__minSolvDensity): + newTube[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] = densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] + if self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] < 0 or self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] == oldGroup: + newVoxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] = oldGroup + else: + tmp = str(self.__voxelGroupGrid[convertPosition[0]][convertPosition[1]][convertPosition[2]+i]) + tmp += " "+str(oldGroup) + newVoxelGroupGrid[convertPosition[0]+i][convertPosition[1]][convertPosition[2]] = tmp + ## adding voxel to groups + #[[x][y][z][x][y][z]] + x = convertPosition[0] + y = convertPosition[1] + z = convertPosition[2]+i + self.__groups[oldGroup].append([x,y,z]) + i -=1 + if (densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] < self.__minSolvDensity): + i = 1 + if (oldTube[convertPosition[0]][convertPosition[1]][convertPosition[2]+i] > 0): + i = 1 + if (densityProtein[xCount][yCount][zCount+i] > float(self.__minProteinDensity)): + i = 1 + for xCount in range (len (newTube)): + for yCount in range (len (newTube[0])): + for zCount in range (len (newTube[0][0])): + if newTube[xCount][yCount][zCount] > 0: + oldTube[xCount][yCount][zCount] = newTube[xCount][yCount][zCount] + self.__voxelGroupGrid[xCount][yCount][zCount] = newVoxelGroupGrid[xCount][yCount][zCount] + self.__newTubeDensity = oldTube diff --git a/lib/filltube.pyc b/lib/filltube.pyc new file mode 100644 index 0000000..7695625 Binary files /dev/null and b/lib/filltube.pyc differ diff --git a/lib/filterdx.py b/lib/filterdx.py new file mode 100644 index 0000000..4072c4b --- /dev/null +++ b/lib/filterdx.py @@ -0,0 +1,167 @@ +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +### +# +# v0.12 +# +FilterByDensity +# +#version 0.066 +#- filter by neighbourmin varriable is now adjustable from gui +#version 0.062 +# changelog +#- protein & water dxObject are not nessasary any more +### +from buildopendx import * +class FilterDx(): + __stepSize = () + __tube = [] + __minDiameter = () + __dxObjectTube =() + + def __init__(self, dxObjectTube, minDiameter=0): + self.__tube = dxObjectTube.getDensity() + self.__dxObjectTube = dxObjectTube + self.__stepSize = dxObjectTube.getStepsize() + self.__minDiameter = int (float(minDiameter) / float(self.__stepSize)) + + def filterTunnelNeighbour(self = None,neighbourMin=10): + print "Start filtering cavities by neighbour "+str(neighbourMin) + tube = [[[0 for zR in range(len(self.__tube[0][0]))] for yR in range(len(self.__tube[0]))] for xR in range(len(self.__tube))] + for xCount in range (len (self.__tube)): + for yCount in range (len (self.__tube[0])): + for zCount in range (len (self.__tube[0][0])): +# counting neighbour voxels with densities + counter = 0 + if (self.__tube[xCount][yCount][zCount] < 0.01): + continue + list = [-1,0,1] + for xPos in list: + for yPos in list: + for zPos in list: + if (self.__tube[xCount+xPos][yCount+yPos][zCount+zPos] > 0.01): + counter +=1 + if (counter > neighbourMin): + tube[xCount][yCount][zCount] = self.__tube[xCount][yCount][zCount] +# outputLine = int ((xCount*100)/(len(self.__tube))) +# print str(outputLine)+" %\b\b\b\b\b", + self.__tube = tube + + def filterDistance(self = None): + print "Start filtering by Distance "+str(self.__minDiameter) + tube = [[[0 for zR in range(len(self.__tube[0][0]))] for yR in range(len(self.__tube[0]))] for xR in range(len(self.__tube))] + for xCount in range (len (self.__tube)): + for yCount in range (len (self.__tube[0])): + for zCount in range (len (self.__tube[0][0])): +# scanning x direction + xDiameter = 1 + for xScan in range (xCount,len (self.__tube)): + if (self.__tube[xScan][yCount][zCount] > 0.01 ): + xDiameter += 1 + if (self.__tube[xScan][yCount][zCount] < 0.01 ): + break + + xScan = xCount + while (xScan > 0): + xScan -=1 + if (self.__tube[xScan][yCount][zCount] > 0.01 ): + xDiameter +=1 + if (self.__tube[xScan][yCount][zCount] < 0.01 ): + break +#scanning y direction + yDiameter = 1 + for yScan in range (yCount,len (self.__tube[0])): + if (self.__tube[xCount][yScan][zCount] > 0.01 ): + yDiameter += 1 + if (self.__tube[xCount][yScan][zCount] < 0.01 ): + break + yScan = yCount + while (yScan > 0): + yScan -=1 + if (self.__tube[xCount][yScan][zCount] > 0.01 ): + yDiameter +=1 + if (self.__tube[xCount][yScan][zCount] < 0.01 ): + break +#scanning z direction + zDiameter = 1 + for zScan in range (zCount,len (self.__tube[0][0])): + if (self.__tube[xCount][yCount][zScan] > 0.01 ): + zDiameter += 1 + if (self.__tube[xCount][yCount][zScan] < 0.01 ): + break + zScan = zCount + while (zScan > 0): + zScan -=1 + if (self.__tube[xCount][yCount][zScan] > 0.01 ): + zDiameter +=1 + if (self.__tube[xCount][yCount][zScan] < 0.01 ): + break + if (xDiameter +yDiameter+zDiameter < 2*self.__minDiameter): + continue +# 2D could be implemented in future +# if (xDiameter > self.__minDiameter and yDiameter > self.__minDiameter and zDiameter > 1): +# tube[xCount][yCount][zCount] = self.__tube[xCount][yCount][zCount] +# if (xDiameter > self.__minDiameter and zDiameter > self.__minDiameter and yDiameter > 1): +# tube[xCount][yCount][zCount] = self.__tube[xCount][yCount][zCount] +# +# if (yDiameter > self.__minDiameter and zDiameter > self.__minDiameter and xDiameter > 1): +# tube[xCount][yCount][zCount] = self.__tube[xCount][yCount][zCount] +# + +#3D + if (xDiameter > self.__minDiameter and yDiameter > self.__minDiameter and zDiameter > self.__minDiameter): + tube[xCount][yCount][zCount] = self.__tube[xCount][yCount][zCount] + + +# outputLine = int ((xCount*100)/(len(self.__tube))) +# print str(outputLine)+" %\b\b\b\b\b", + self.__tube = tube + + def filterByDensity(self,densityThreshold): + print "Minimum residence probability: "+densityThreshold + print "Start filtering..." + densities = self.__dxObjectTube.getDensity() + densities_filtered = [[[-1 for zR in range(len (densities[0][0]))] for yR in range(len (densities[0]))] for xR in range(len (densities))] + for xCount in range (len (densities)): + for yCount in range (len (densities[0])): + for zCount in range (len (densities[0][0])): + if float(densities[xCount][yCount][zCount]) > float(densityThreshold): + densities_filtered[xCount][yCount][zCount] = densities[xCount][yCount][zCount] + outputLine = int ((xCount*100)/len(densities)) + print str(outputLine)+" %\r", +# tmp = (BuildOpenDx (density = densities_filtered, stepSize = tubeDx.getDxObject().getStepsize(), origin = tubeDx.getDxObject().getOrigin(),dimention = tubeDx.getDxObject().getDimention())) + self.__dxObjectTube.setDensity(densities_filtered) + +# This routine deletes all densities INSIDE the given rectangle defined by (xmin,ymin,zmin) and (xmax,ymax,zmax) UNTESTET ROUTINE BEWARE !!!!!!! + def deleteRectangleInside(self,minArray,mayArray): + densities = self.__dxObjectTube.getDensity() +# densities_filtered = [[[-1 for zR in range(len (densities[0][0]))] for yR in range(len (densities[0]))] for xR in range(len (densities))] + origin = self.__dxObjectTube.getOrigin() + for xCount in range (len (densities)): + for yCount in range (len (densities[0])): + for zCount in range (len (densities[0][0])): + x = float (float(xCount) * float(self.__stepSize)) + float(origin[0]) + y = float (float(yCount) * float(self.__stepSize)) + float(origin[1]) + z = float (float(zCount) * float(self.__stepSize)) + float(origin[2]) + if x > float(minArray[0]) and x < float(maxArray[0]): + if y > float(minArray[1]) and y < float(maxArray[1]): + if z > float(minArray[2]) and z < float(maxArray[2]): + densities[xCount][yCount][zCount] = 0 + self.__dxObjectTube.setDensity(densities) + + def getDxObject(self=None): + object = BuildOpenDx (self.__dxObjectTube.getDensity(), self.__dxObjectTube.getStepsize(), self.__dxObjectTube.getOrigin(), self.__dxObjectTube.getDimention()) + return object diff --git a/lib/filterdx.pyc b/lib/filterdx.pyc new file mode 100644 index 0000000..735fa54 Binary files /dev/null and b/lib/filterdx.pyc differ diff --git a/lib/grouptubes.py b/lib/grouptubes.py new file mode 100644 index 0000000..dd3bc47 --- /dev/null +++ b/lib/grouptubes.py @@ -0,0 +1,208 @@ +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +##### +# v0.15 +# -excluded diagonals from grouping +# + + +import sys +class GroupTubes(): + + __tubeDensity = [] + __tubeGroup =[] +# tubeVoxelGrouped =[] + __groupCount = 0 + __groupVoxel= [] + def __init__(self, dxObjectTube): + self.__tubeDensity = dxObjectTube.getDensity() + + def testIndex(self, string, array): + try: + i = array.index(string) + except ValueError: + i = -1 # string not found + return i + + def findGroups(self): + print "Start grouping cavities" +#in tubeVoxelGrouped is [x][y][z] = group + tubeVoxelGrouped =[[[-1 for zR in range(len (self.__tubeDensity[0][0]))] for yR in range(len (self.__tubeDensity[0]))] for xR in range(len (self.__tubeDensity))] +#self.__tubeGroup : +#[[x][y][z][x][y][z]] + self.__tubeGroup = [] + for xCount in range (len (self.__tubeDensity)): + for yCount in range (len (self.__tubeDensity[0])): + for zCount in range (len (self.__tubeDensity[0][0])): +#analysing neighbourhood if there is allready a cluster: + if (self.__tubeDensity[xCount][yCount][zCount]): + list = [-1,0] + zlist = [-1,0,1] + foundGroup = -1 + for xPos in list: + for yPos in list: + for zPos in zlist: + if (zPos == 0 and xPos == 0 and yPos == 0): + continue +#not visited yet + if (xPos == 0 and yPos == 0 and zPos == 1): + continue + +### excluding diagonals EXPERIMENTAL ############## +# if (xPos == -1 and yPos == -1 and zPos == -1): +# continue +# if (xPos == 0 and yPos == -1 and zPos == -1): +# continue +# if (xPos == -1 and yPos == -1 and zPos == 0 ): +# continue +# if (xPos == -1 and yPos == -1 and zPos == +1 ): +# continue +# if (xPos == -1 and yPos == 0 and zPos == +1): +# continue +# if (xPos == -1 and yPos == 0 and zPos == -1 ): +# continue +######################################### + + +#test if another group is defined in neighbourhood: + if (self.__tubeDensity[xCount+xPos][yCount+yPos][zCount+zPos] ): + coordinate = (xCount+xPos,yCount+yPos,zCount+zPos) + neighbourGrp = tubeVoxelGrouped [coordinate[0]][coordinate[1]][coordinate[2]] + + if (foundGroup >= 0): +# neighbourGrp = tubeVoxelGrouped [coordinate[0]][coordinate[1]][coordinate[2]] + if (neighbourGrp == foundGroup): + continue +# if another group is found, members will be copied from later found cluster to earlier found cluster +# print "merging two groups: "+ str(neighbourGrp) + " " + str(foundGroup) + if (neighbourGrp > foundGroup): + + for x in range (len (tubeVoxelGrouped)): + for y in range (len (tubeVoxelGrouped[0])): + for z in range (len (tubeVoxelGrouped[0][0])): + if (tubeVoxelGrouped [x][y][z] > neighbourGrp): + tubeVoxelGrouped [x][y][z] -= 1 +#test if a goup has only one member +# member member member +# [[x][y][z][x][y][z]] OR [[x][y][z]] +# G R O U P G R O U P + +# a try catch block is nessasary to decide if a Group only has one member +# if a group has more than one member +# member gets his own dimension +# if a group is only build up by one member is has not the dimension of member + + try: + x = self.__tubeGroup[neighbourGrp][0][0] + oneMember = 0 + except TypeError: + x = self.__tubeGroup[neighbourGrp][0] + y = self.__tubeGroup[neighbourGrp][1] + z = self.__tubeGroup[neighbourGrp][2] + self.__tubeGroup[foundGroup].append([x,y,z]) + oneMember = 1 + tubeVoxelGrouped [x][y][z] = foundGroup + + if (oneMember == 0): + for member in range (len(self.__tubeGroup[neighbourGrp])): + x = self.__tubeGroup[neighbourGrp][member][0] + y = self.__tubeGroup[neighbourGrp][member][1] + z = self.__tubeGroup[neighbourGrp][member][2] + tubeVoxelGrouped [x][y][z] = foundGroup + self.__tubeGroup[foundGroup].extend(self.__tubeGroup[neighbourGrp]) + + + del self. __tubeGroup[neighbourGrp] + neighbourGrp = foundGroup + continue + if (neighbourGrp < foundGroup): + for x in range (len (tubeVoxelGrouped)): + for y in range (len (tubeVoxelGrouped[0])): + for z in range (len (tubeVoxelGrouped[0][0])): + if (tubeVoxelGrouped [x][y][z] > foundGroup): + tubeVoxelGrouped [x][y][z] -= 1 + try: + x = self.__tubeGroup[foundGroup][0][0] + oneMember = 0 + except TypeError: + x = self.__tubeGroup[foundGroup][0] + y = self.__tubeGroup[foundGroup][1] + z = self.__tubeGroup[foundGroup][2] + self.__tubeGroup[neighbourGrp].append([x,y,z]) + oneMember = 1 + tubeVoxelGrouped [x][y][z] = neighbourGrp + if (oneMember == 0): + for member in range (len(self.__tubeGroup[foundGroup])): + x = self.__tubeGroup[foundGroup][member][0] + y = self.__tubeGroup[foundGroup][member][1] + z = self.__tubeGroup[foundGroup][member][2] + tubeVoxelGrouped [x][y][z] = neighbourGrp + self.__tubeGroup[neighbourGrp].extend(self.__tubeGroup[foundGroup]) + + del self.__tubeGroup[foundGroup] + foundGroup = neighbourGrp + continue +# searching for an existing group in neighbourhood + if (foundGroup < 0): + self.__tubeGroup[neighbourGrp].append([xCount,yCount,zCount]) + foundGroup = neighbourGrp + tubeVoxelGrouped [xCount][yCount][zCount] = foundGroup + + +# if no group was found in neighbourhood a new group will be created + if (foundGroup == -1 ): +# print "create new group: "+ str(len(self.__tubeGroup)) + foundGroup = len(self.__tubeGroup) + tupelTMP = [([xCount,yCount,zCount])] + self.__tubeGroup.append (tupelTMP) + tubeVoxelGrouped [xCount][yCount][zCount] = foundGroup + +# outputLine = int ((xCount*100)/(len(self.__tubeDensity))) + outputLine = int ((xCount*100)/(len(self.__tubeDensity))) +# print str(outputLine)+" %\b\b\b\b\b", + print str(outputLine)+" %\r", + sys.stdout.flush() + print "final groups: "+str(len(self.__tubeGroup)) + self.__groupVoxel = tubeVoxelGrouped +#returns +# member member +# [[x][y][z][x][y][z]] +# G R O U P +# groupes in TubeDx + def getTubeGroup(self): + return self.__tubeGroup + +#groupes +#[[x][y][z][x][y][z]] + + def getGroupes (self): + return self.__tubeGroup +#### =====:-( whoa scary those getters are + +#voxelGroupGrid in tubeDX +# returns [x][y][z] = group + def getGroups(self): + return self.__groupVoxel + + + +#voxelGroupGrid +#[x][y][z] = group + + def getVoxelGroupGrid (self): + return self.__groupVoxel + diff --git a/lib/grouptubes.pyc b/lib/grouptubes.pyc new file mode 100644 index 0000000..f8a20bf Binary files /dev/null and b/lib/grouptubes.pyc differ diff --git a/lib/planesearch.py b/lib/planesearch.py new file mode 100644 index 0000000..e332248 --- /dev/null +++ b/lib/planesearch.py @@ -0,0 +1,469 @@ +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + +from buildopendx import * +##### +# +# +# v0.2 +# - __storeOverLappingISV = False <--- switches if ISV will be deleted if a protein Voxels > protMinDens is located at the very same position +# +# +# v0.15 +# +new method " +# -minSolventDens(default)= 0 now +# +# v.0.091 +# +support solvent files smaller than protein files +# +# v.0.08 +# +user adjustable solvent density + + +class PlaneSearch(): + __densityWater = [] + __densityProtein = [] + __stepSize = () + __tube = [] + __dxObjectWater =() + __dxObjectProtein =() + __minDiameter = () + __minProteinDensity =() + __secondDensities= False +### Trigger that switches between old mode ( deleting ISV overlapping with protein Voxels) and new mode (ignore overlap) + __storeOverLappingISV = False + + def __init__(self, dxObjectWater, dxObjectProtein, minDiameter="0",minProteinDens = "0.005",minSolventDens = "0.00", mdk = False, storeOverLappingISV = False): + self.__minProteinDensity = float(minProteinDens) + self.__densityWater = dxObjectWater.getDensity() + self.__stepSize = dxObjectWater.getStepsize() + self.__minDiameter = int (float(minDiameter) / float(self.__stepSize)) + print 'Protein Threshold = ' + minProteinDens+'\nSolvent Threshold = ' + minSolventDens + print 'minDiameter: '+ str(minDiameter)+ ' Angstrom \nVolmap Stepsize: '+str(self.__stepSize)+ '\nResulting minimum voxel minDiameter: ' +str(self.__minDiameter) +# self.__minDiameter = minDiameter +#implementing solvent user adjustable + self.__minSolventDens= float(minSolventDens) + self.__densityProtein = dxObjectProtein.getDensity() + self.__dxObjectWater = dxObjectWater + self.__dxObjectProtein = dxObjectProtein + self.__tube = [[[0 for zR in range(len (self.__densityWater[0][0]))] for yR in range(len (self.__densityWater[0]))] for xR in range(len (self.__densityWater))] + self.__mdk = mdk + self.__storeOverLappingISV = storeOverLappingISV +# plane must be 'x' 'y' or 'z' + def oneDimension(self, plane): + proteinOrigin = self.__dxObjectProtein.getOrigin() +#stepsize multiplier need to be implemented if stepsize !=1 + stepMultiplier = 1 / self.__stepSize + waterOrigin = self.__dxObjectWater.getOrigin() + waterToProteinOrigin = (int ((proteinOrigin[0] - waterOrigin[0])*stepMultiplier),int ((proteinOrigin[1] - waterOrigin[1])*stepMultiplier),int ((proteinOrigin[2] - waterOrigin[2])*stepMultiplier)) + print 'Scanning '+ plane +' plane minimum diameter is set to: '+ str (self.__minDiameter)+ ' voxels' + row = () + if (plane == 'z'): + for xCount in range (len (self.__densityProtein)): + for yCount in range (len (self.__densityProtein[0])): + row =[] + for zCount in range (len (self.__densityProtein[0][0])): + row.append (self.__densityProtein[xCount][yCount][zCount]) + for i in range (len(row)): + if (self.__densityProtein[xCount][yCount][i] > self.__minProteinDensity and not self.__storeOverLappingISV ): + continue + convertPosition = (xCount+waterToProteinOrigin[0], yCount+waterToProteinOrigin[1],i+waterToProteinOrigin[2]) + try: + self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + except IndexError: + continue +## implementing solvent variables: + if (self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] < self.__minSolventDens): +# if (self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] < 0.01): + continue + proteinLeft = 0 + proteinRight = 0 + + j= i + + while (j > 0): + j -=1 + if (self.__densityProtein[xCount][yCount][j] >= self.__minProteinDensity): + proteinLeft = i - j + break + + for k in range (i,len (row)): + if (self.__densityProtein[xCount][yCount][k] >= self.__minProteinDensity ): + proteinRight =k-i + break + borderSum = proteinLeft+proteinRight +# print 'borderSum: ' + str(borderSum) + if (proteinLeft > 0 and proteinRight > 0 and borderSum > self.__minDiameter ): + self.__tube[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + + + if (plane == 'y'): + for xCount in range (len (self.__densityProtein)): + for zCount in range (len (self.__densityProtein[0][0])): + row =[] + for yCount in range (len (self.__densityProtein[0])): + row.append (self.__densityProtein[xCount][yCount][zCount]) + for i in range (len(row)): + if (self.__densityProtein[xCount][i][zCount] > self.__minProteinDensity and not self.__storeOverLappingISV): + continue + convertPosition = (xCount+waterToProteinOrigin[0], i+waterToProteinOrigin[1],zCount+waterToProteinOrigin[2]) + try: + self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + except IndexError: + continue +## implementing solvent variables: + if (self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] < self.__minSolventDens): +# if (self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] < 0.01): + continue +# testing if water lies inside the protein + proteinLeft = 0 + proteinRight = 0 + j= i + while (j > 0): + j -=1 + if (self.__densityProtein[xCount][j][zCount] >= self.__minProteinDensity ): + proteinLeft = i - j + break + + for k in range (i,len (row)): + if (self.__densityProtein[xCount][k][zCount] >= self.__minProteinDensity ): + proteinRight = k - i + break + borderSum = proteinLeft+proteinRight + if (proteinLeft > 0 and proteinRight > 0 and borderSum > self.__minDiameter ): + self.__tube[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + + if (plane == 'x'): + for yCount in range (len (self.__densityProtein[0])): + for zCount in range (len (self.__densityProtein[0][0])): + row =[] + for xCount in range (len (self.__densityProtein)): + row.append (self.__densityProtein[xCount][yCount][zCount]) + for i in range (len(row)): + if (self.__densityProtein[i][yCount][zCount] > self.__minProteinDensity and not self.__storeOverLappingISV): + continue + convertPosition = (i+waterToProteinOrigin[0], yCount+waterToProteinOrigin[1],zCount+waterToProteinOrigin[2]) + try: + self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + except IndexError: + continue +## implementing solvent variables: + if (self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] < self.__minSolventDens): +# if (self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] < 0.01): + continue +# testing if water lies inside the protein + proteinLeft = 0 + proteinRight = 0 + j= i + while (j > 0): + j -=1 + if (self.__densityProtein[j][yCount][zCount] >= self.__minProteinDensity ): + proteinLeft = i - j + break + + for k in range (i,len (row)): + if (self.__densityProtein[k][yCount][zCount] >= self.__minProteinDensity ): + proteinRight = k-i + break + borderSum = proteinLeft+proteinRight + if (proteinLeft > 0 and proteinRight > 0 and borderSum > self.__minDiameter): + self.__tube[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + print "Done" + + def twoDimension(self = None): + print "Scanning 2D" + proteinOrigin = self.__dxObjectProtein.getOrigin() + waterOrigin = self.__dxObjectWater.getOrigin() +#stepsize multiplier need to be implemented if stepsize !=1 + stepMultiplier = 1 / self.__stepSize + waterToProteinOrigin = (int ((proteinOrigin[0] - waterOrigin[0])*stepMultiplier),int ((proteinOrigin[1] - waterOrigin[1])*stepMultiplier),int ((proteinOrigin[2] - waterOrigin[2])*stepMultiplier)) + for xCount in range (len (self.__densityProtein)): + for yCount in range (len (self.__densityProtein[0])): + for zCount in range (len (self.__densityProtein[0][0])): + + convertPosition = (xCount+waterToProteinOrigin[0], yCount+waterToProteinOrigin[1],zCount+waterToProteinOrigin[2]) + try: + self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + except IndexError: + continue + if (self.__densityProtein[xCount][yCount][zCount] > self.__minProteinDensity and not self.__storeOverLappingISV): + continue + xProteinLeft = 0 + xProteinRight = 0 + j = xCount + while (j > 0): + j -=1 + if (self.__densityProtein[j][yCount][zCount] >= self.__minProteinDensity): + xProteinLeft = xCount - j + break + + for k in range (xCount+1,len(self.__densityProtein)): #### xcount +1 !!! to test neighbours and not the voxel itself !!! + if (self.__densityProtein[k][yCount][zCount] >= self.__minProteinDensity ): + xProteinRight = k- xCount + break + xBorderSum = xProteinLeft+xProteinRight + + yProteinLeft = 0 + yProteinRight = 0 + j = yCount + while (j > 0): + j -=1 + if (self.__densityProtein[xCount][j][zCount] >= self.__minProteinDensity ): + yProteinLeft = yCount - j + break + + for k in range (yCount+1,len(self.__densityProtein[0])): + if (self.__densityProtein[xCount][k][zCount] >= self.__minProteinDensity ): + yProteinRight = k - yCount + break + yBorderSum = yProteinLeft+yProteinRight + + zProteinLeft = 0 + zProteinRight = 0 + j = zCount + while (j > 0): + j -=1 + if (self.__densityProtein[xCount][yCount][j] >= self.__minProteinDensity ): + zProteinLeft = zCount - j + break + + for k in range (zCount+1,len(self.__densityProtein[0][0])): + if (self.__densityProtein[xCount][yCount][k] >= self.__minProteinDensity ): + zProteinRight =k- zCount + break + zBorderSum = zProteinLeft+zProteinRight + if (xProteinLeft > 0 and xProteinRight > 0 ): + if (yProteinLeft > 0 and yProteinRight > 0): + if (xBorderSum > self.__minDiameter and yBorderSum > self.__minDiameter and self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] > self.__minSolventDens): +# if (xBorderSum > self.__minDiameter and yBorderSum > self.__minDiameter): + self.__tube[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + if (zProteinLeft > 0 and zProteinRight > 0): + if (xBorderSum > self.__minDiameter and zBorderSum > self.__minDiameter and self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] > self.__minSolventDens): + self.__tube[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + + if (yProteinLeft > 0 and yProteinRight > 0 and zProteinLeft > 0 and zProteinRight > 0 ): + if (yBorderSum > self.__minDiameter and zBorderSum > self.__minDiameter and self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] > self.__minSolventDens): +# if (yBorderSum > self.__minDiameter and zBorderSum > self.__minDiameter): + self.__tube[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + outputLine = int ((xCount*100)/len(self.__densityProtein)) + print str(outputLine)+" %\b\b\b\b\b", + print "Done" +## +## This methods saves isv below a required minDiameter neighvborhood in a separate density file +## 2 dxobjects will be created and can be procecced via grouping separatly +## + def twoDimensionMinDiaArray(self = None): + print "Scanning 2D \nKeeping ISV < " + str(self.__minDiameter)+ " Angstrom" + proteinOrigin = self.__dxObjectProtein.getOrigin() + waterOrigin = self.__dxObjectWater.getOrigin() + + belowMinDiameterArray = [[[0 for zR in range(len (self.__densityWater[0][0]))] for yR in range(len (self.__densityWater[0]))] for xR in range(len (self.__densityWater))] +#stepsize multiplier need to be implemented if stepsize !=1 + stepMultiplier = 1 / self.__stepSize + waterToProteinOrigin = (int ((proteinOrigin[0] - waterOrigin[0])*stepMultiplier),int ((proteinOrigin[1] - waterOrigin[1])*stepMultiplier),int ((proteinOrigin[2] - waterOrigin[2])*stepMultiplier)) + for xCount in range (len (self.__densityProtein)): + for yCount in range (len (self.__densityProtein[0])): + for zCount in range (len (self.__densityProtein[0][0])): + + convertPosition = (xCount+waterToProteinOrigin[0], yCount+waterToProteinOrigin[1],zCount+waterToProteinOrigin[2]) + try: + self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + except IndexError: + continue + if (self.__densityProtein[xCount][yCount][zCount] > self.__minProteinDensity and not self.__storeOverLappingISV): + continue + xProteinLeft = 0 + xProteinRight = 0 + j = xCount + while (j > 0): + j -=1 + if (self.__densityProtein[j][yCount][zCount] >= self.__minProteinDensity): + xProteinLeft = xCount - j + break + + for k in range (xCount+1,len(self.__densityProtein)): + if (self.__densityProtein[k][yCount][zCount] >= self.__minProteinDensity ): + xProteinRight = k- xCount + break + xBorderSum = xProteinLeft+xProteinRight + + yProteinLeft = 0 + yProteinRight = 0 + j = yCount + while (j > 0): + j -=1 + if (self.__densityProtein[xCount][j][zCount] >= self.__minProteinDensity ): + yProteinLeft = yCount - j + break + + for k in range (yCount+1,len(self.__densityProtein[0])): + if (self.__densityProtein[xCount][k][zCount] >= self.__minProteinDensity ): + yProteinRight = k - yCount + break + yBorderSum = yProteinLeft+yProteinRight + + zProteinLeft = 0 + zProteinRight = 0 + j = zCount + while (j > 0): + j -=1 + if (self.__densityProtein[xCount][yCount][j] >= self.__minProteinDensity ): + zProteinLeft = zCount - j + break + + for k in range (zCount+1,len(self.__densityProtein[0][0])): + if (self.__densityProtein[xCount][yCount][k] >= self.__minProteinDensity ): + zProteinRight =k- zCount + break + zBorderSum = zProteinLeft+zProteinRight + + if (xProteinLeft > 0 and xProteinRight > 0 and yProteinLeft > 0 and yProteinRight > 0): + if (xBorderSum >= self.__minDiameter and yBorderSum >= self.__minDiameter and self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] > self.__minSolventDens): + self.__tube[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] +# +# keeping isv below required diameter in a separete Array +# v0.2 + elif self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] > self.__minSolventDens: + belowMinDiameterArray[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] +############################################# + if (xProteinLeft > 0 and xProteinRight > 0 and zProteinLeft > 0 and zProteinRight > 0): + if (xBorderSum >= self.__minDiameter and zBorderSum >= self.__minDiameter and self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] > self.__minSolventDens): + self.__tube[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + elif self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] > self.__minSolventDens: + belowMinDiameterArray[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + + if (yProteinLeft > 0 and yProteinRight > 0 and zProteinLeft > 0 and zProteinRight > 0 ): + if (yBorderSum >= self.__minDiameter and zBorderSum >= self.__minDiameter and self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] > self.__minSolventDens): + self.__tube[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + elif self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] > self.__minSolventDens: + belowMinDiameterArray[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + + +### +# mdk block +### +# if (xProteinLeft > 0 and xProteinRight > 0) and (yProteinLeft > 0 and yProteinRight > 0): +# if ((xBorderSum <= self.__minDiameter or yBorderSum <= self.__minDiameter)# or (xBorderSum < self.__minDiameter or yBorderSum <= self.__minDiameter) +# and self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] > self.__minSolventDens): +# belowMinDiameterArray[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] +# +## if (xProteinLeft > 0 and xProteinRight > 0) and (zProteinLeft > 0 and zProteinRight > 0): +# if ((xBorderSum <= self.__minDiameter or zBorderSum <= self.__minDiameter) #or (xBorderSum < self.__minDiameter or zBorderSum <= self.__minDiameter) +# +# and self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] > self.__minSolventDens): +# belowMinDiameterArray[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[c#onvertPosition[0]][convertPosition[1]][convertPosition[2]] +# +# if (yProteinLeft > 0 and yProteinRight > 0) and (zProteinLeft > 0 and zProteinRight > 0): +# if ((yBorderSum <= self.__minDiameter or zBorderSum <= self.__minDiameter)# or (yBorderSum < self.__minDiameter or zBorderSum <= self.__minDiameter) +# and self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] > self.__minSolventDens): +# belowMinDiameterArray[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + + outputLine = int ((xCount*100)/len(self.__densityProtein)) + print str(outputLine)+" %\b\b\b\b\b", + self.__secondDensities = belowMinDiameterArray + print "Done" + + + + def threeDimension(self = None): + print "Scanning 3D" +#stepsize multiplier need to be implemented if stepsize !=1 + stepMultiplier = 1 / self.__stepSize + proteinOrigin = self.__dxObjectProtein.getOrigin() + waterOrigin = self.__dxObjectWater.getOrigin() + + if self.__mdk: + belowMinDiameterArray = [[[0 for zR in range(len (self.__densityWater[0][0]))] for yR in range(len (self.__densityWater[0]))] for xR in range(len (self.__densityWater))] + + waterToProteinOrigin = (int ((proteinOrigin[0] - waterOrigin[0])*stepMultiplier),int ((proteinOrigin[1] - waterOrigin[1])*stepMultiplier),int ((proteinOrigin[2] - waterOrigin[2])*stepMultiplier)) + for xCount in range (len (self.__densityProtein)): + for yCount in range (len (self.__densityProtein[0])): + for zCount in range (len (self.__densityProtein[0][0])): + convertPosition = (xCount+waterToProteinOrigin[0], yCount+waterToProteinOrigin[1],zCount+waterToProteinOrigin[2]) + try: + self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + except IndexError: + continue + + if (self.__densityProtein[xCount][yCount][zCount] > self.__minProteinDensity and not self.__storeOverLappingISV): + continue + xProteinLeft = 0 + xProteinRight = 0 + j = xCount + while (j > 0): + j -=1 + if (self.__densityProtein[j][yCount][zCount] >= self.__minProteinDensity): + xProteinLeft = xCount - j + break + + for k in range (xCount+1,len(self.__densityProtein)): + if (self.__densityProtein[k][yCount][zCount] >= self.__minProteinDensity ): + xProteinRight = k- xCount + break + xBorderSum = xProteinLeft+xProteinRight + + yProteinLeft = 0 + yProteinRight = 0 + j = yCount + while (j > 0): + j -=1 + if (self.__densityProtein[xCount][j][zCount] >= self.__minProteinDensity ): + yProteinLeft = yCount - j + break + + for k in range (yCount+1,len(self.__densityProtein[0])): + if (self.__densityProtein[xCount][k][zCount] >= self.__minProteinDensity ): + yProteinRight = k - yCount + break + yBorderSum = yProteinLeft+yProteinRight + + zProteinLeft = 0 + zProteinRight = 0 + j = zCount + while (j > 0): + j -=1 + if (self.__densityProtein[xCount][yCount][j] >= self.__minProteinDensity ): + zProteinLeft = zCount - j + break + + for k in range (zCount+1,len(self.__densityProtein[0][0])): + if (self.__densityProtein[xCount][yCount][k] >= self.__minProteinDensity ): + zProteinRight =k- zCount + break + zBorderSum = zProteinLeft+zProteinRight + if (xProteinLeft > 0 and xProteinRight > 0 and yProteinLeft > 0 and yProteinRight > 0 and zProteinLeft > 0 and zProteinRight > 0 ): + if (xBorderSum >= self.__minDiameter and yBorderSum >= self.__minDiameter and zBorderSum >= self.__minDiameter and self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] > self.__minSolventDens): +# if (xBorderSum > self.__minDiameter and xBorderSum > self.__minDiameter and zBorderSum > self.__minDiameter): + self.__tube[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + elif self.__mdk: + if self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] > self.__minSolventDens: + belowMinDiameterArray[convertPosition[0]][convertPosition[1]][convertPosition[2]] = self.__densityWater[convertPosition[0]][convertPosition[1]][convertPosition[2]] + outputLine = int ((xCount*100)/len(self.__densityProtein)) + print str(outputLine)+" %\b\b\b\b\b", + + if self.__mdk: + self.__secondDensities = belowMinDiameterArray + print "Done" + def getDxObject(self=None): + object = BuildOpenDx (self.__tube, self.__dxObjectWater.getStepsize(), self.__dxObjectWater.getOrigin(),self.__dxObjectWater.getDimention()) + return object + def getDxObjectArray(self=None): + object1 = BuildOpenDx (self.__tube, self.__dxObjectWater.getStepsize(), self.__dxObjectWater.getOrigin(),self.__dxObjectWater.getDimention()) + object2 = BuildOpenDx (self.__secondDensities, self.__dxObjectWater.getStepsize(), self.__dxObjectWater.getOrigin(),self.__dxObjectWater.getDimention()) + object = [object1,object2] + return object + diff --git a/lib/planesearch.pyc b/lib/planesearch.pyc new file mode 100644 index 0000000..a41abd5 Binary files /dev/null and b/lib/planesearch.pyc differ diff --git a/lib/postgroup.py b/lib/postgroup.py new file mode 100644 index 0000000..a1c7c7f --- /dev/null +++ b/lib/postgroup.py @@ -0,0 +1,665 @@ + +# +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +# +# +# v0.23 +# +'core_cavities_after_filter.pdb' stores the grouped core cavities +# +# dxTuber v0.2 +# +new routine to postgroup / seperate existing cavities, deviding cavities in subcavities +# +# + +import re, math +from tubedx import * +from grouptubes import * +from buildopendx import * +from settings import * +### for testing only ... +from writepdb import * +#### + + + + + +#### this class seperates previously detected cavities. Postgrouping. +class PostGroup(): +#groupes +#[[x][y][z][x][y][z]] + __cavityGroup =[] +# tubeVoxelGrouped =[] + __groupCount = 0 + + +# returns [x][y][z] = group + __groupVoxel= [] + +## opendx Object containing the cavity voxels + __cavitysDxObject = () + +# the input tubedx + __tubeDx = () + __version = () + def __init__(self, tubeDx): +# returns [x][y][z] = group + self.__groupVoxel = tubeDx.getVoxelGroupGrid() +#groupes +#[[x][y][z][x][y][z]] + self.__cavityGroup = tubeDx.getGroupes() + + self.__cavitysDxObject = tubeDx.getDxObject() + self.__tubeDx = tubeDx + settings = Settings() + self.__version = settings.getVersion() + + +## postgroupArray contains a collection of cavityID's to postgroup [int] +#returns an "density" array containung Neighborhod filtered values + def calcNeighbors (self, cavityID, edgeLength = 2 ): + +# cubicSize = (edgeLength+edgeLength+1 ) * (edgeLength+edgeLength+1 ) * (edgeLength+edgeLength+1 ) #fix this !!! +# print "Calculating for each voxel the sum of its neighbors in a "+ str(cubicSize) +" cubic box" # fix this ... + + boxArray = [[[-1 for zR in range(len (self.__cavitysDxObject.getDensity()[0][0]))] for yR in range(len (self.__cavitysDxObject.getDensity()[0]))] for xR in range(len (self.__cavitysDxObject.getDensity()))] +# for cavityID in postGroupArray: + for x in range (len (self.__cavitysDxObject.getDensity())): + for y in range (len (self.__cavitysDxObject.getDensity()[0])): + for z in range (len (self.__cavitysDxObject.getDensity()[0][0])): +# if not re.match(cavityID , str(self.__groupVoxel[x][y][z])) : +# if not ( cavityID == str(self.__groupVoxel[x][y][z]) ): + +## if not cavity : all cavities will be regrouped :) + if cavityID != 'all': + if not ( cavityID == str(self.__groupVoxel[x][y][z]) ): + continue + currentBox = 1 + negEdgeLenght = edgeLength - edgeLength - edgeLength + for xCurrent in range (negEdgeLenght, (edgeLength+1)): + for yCurrent in range (negEdgeLenght, (edgeLength+1)): + for zCurrent in range (negEdgeLenght, (edgeLength+1)): + + if xCurrent == 0 and yCurrent == 0 and zCurrent == 0 : + continue + try: + if cavityID == str(self.__groupVoxel[x+xCurrent][y+yCurrent][z+zCurrent]) : + currentBox += 1 + elif cavityID == 'all' and self.__cavitysDxObject.getDensity()[x+xCurrent][y+yCurrent][z+zCurrent] > 0: + currentBox += 1 + except IndexError: +# Donothing + a = 1 +# wuha nothing done :> + boxArray[x][y][z] = currentBox + return boxArray + +## +# devides current dxObject into 2 groups containing voxels below (<) threshold and (>=) greater or equal the threshold + def separateByValue(self, threshold, dxObject = __cavitysDxObject): + firstGroup = [[[0 for zR in range(len (dxObject.getDensity()[0][0]))] for yR in range(len (dxObject.getDensity()[0]))] for xR in range(len (dxObject.getDensity()))] + secondGroup = [[[0 for zR in range(len (dxObject.getDensity()[0][0]))] for yR in range(len (dxObject.getDensity()[0]))] for xR in range(len (dxObject.getDensity()))] + for x in range (len (dxObject.getDensity())): + for y in range (len (dxObject.getDensity()[0])): + for z in range (len (dxObject.getDensity()[0][0])): + if dxObject.getDensity()[x][y][z] < threshold: + secondGroup[x][y][z] = dxObject.getDensity()[x][y][z] + else: + firstGroup[x][y][z] = dxObject.getDensity()[x][y][z] + tmpArray = [firstGroup, secondGroup] + return tmpArray +# return firstGroup, secondGroup + +### +# it seems that this method has an malfunction +### + def joinByDistance (self, firstTubeDx, secondOpenDx, distCutoff = -1): + +# first group second group ..... +#[ [x][y][z] [x][y][z] [x][y][z] ] + newGroupes = firstTubeDx.getGroupes() + +#[x][y][z] = groupID + newVoxelGroupGrid = firstTubeDx.getVoxelGroupGrid() + firstDensities = firstTubeDx.getDxObject().getDensity() + + for x in range (len (self.__cavitysDxObject.getDensity())): + for y in range (len (self.__cavitysDxObject.getDensity()[0])): + for z in range (len (self.__cavitysDxObject.getDensity()[0][0])): + if firstTubeDx.getDxObject().getDensity()[x][y][z] > 0: # voxel is member of the "major" (first) selection + continue + size = 1 + foundOneGroup = False + foundOneGroupXYZ = [] + + foundMultipleGroups = [] + foundMultipleGroupsXYZ = [] +## +# +# scanning neighborhood for major cavities .... +## + + while (size <= distCutoff and foundOneGroup == False): +# if distCutoff > 0 and distCutoff > size: +# +# a routine which gathers voxels obove cutoff range should be implemented +# + negLen = size -size -size + posLen = size + + # searching for "major" cavities until the cutOff is reached + for xCurrent in range (negLen, (posLen+1)): + for yCurrent in range (negLen, (posLen+1)): + for zCurrent in range (negLen, (posLen+1)): + + if (x+xCurrent < len(firstTubeDx.getDxObject().getDensity()) and + y+yCurrent < len (firstTubeDx.getDxObject().getDensity()[0]) and + z+zCurrent < len( firstTubeDx.getDxObject().getDensity()[0][0]) and # checking if voxel is in array ... ... index out of bound error ... + x+xCurrent >= 0 and y+yCurrent >= 0 and z+zCurrent >= 0 and + self.__tubeDx.getDxObject().getDensity()[x+xCurrent][y+yCurrent][z+zCurrent] >= self.__tubeDx.getSolventThreshold() # only voxels above previosly defined threshold !! + + ): +# print "do" + +# if firstTubeDx.getDxObject().getDensity()[x+xCurrent][y+yCurrent][z+zCurrent] > 0 and not foundOneGroup: +# foundOneGroupXYZ = [x+xCurrent, y+yCurrent, z+zCurrent] +# foundOneGroup = firstTubeDx.getVoxelGroupGrid()[x+xCurrent][y+yCurrent][z+zCurrent] +# print "yo" +# continue +# if firstTubeDx.getDxObject().getDensity()[x+xCurrent][y+yCurrent][z+zCurrent] > 0 and foundOneGroup: +# foundMultipleGroupsXYZ.append([x+xCurrent, y+yCurrent, z+zCurrent]) +# foundMultipleGroupsXYZ.append(firstTubeDx.getVoxelGroupGrid()[x+xCurrent][y+yCurrent][z+zCurrent]) +# print "yo2" +# print foundMultipleGroupsXYZ + # if firstTubeDx.getDxObject().getDensity()[x+xCurrent][y+yCurrent][z+zCurrent] > 0: + if firstDensities[x+xCurrent][y+yCurrent][z+zCurrent] > 0: + tmpArray = [x+xCurrent, y+yCurrent, z+zCurrent] + foundMultipleGroupsXYZ.append(tmpArray) + foundMultipleGroups.append(firstTubeDx.getVoxelGroupGrid()[x+xCurrent][y+yCurrent][z+zCurrent]) + + + size +=1 +# print "finished while" +#### Evaluating to which major group the current voxel belongs : + + # if foundOneGroup: +# print "found" + if len(foundMultipleGroups) > 1: +# do something + #a= 1 +# inserting current voxel into existing group: +# compare all groups get the smalest distance vector ! +# +# +# print "Ja" + minLength = -1 + minGroup = -1 + for member in range (len(foundMultipleGroupsXYZ)): + xG = foundMultipleGroupsXYZ[member][0] + yG = foundMultipleGroupsXYZ[member][1] + zG = foundMultipleGroupsXYZ[member][2] + length = math.sqrt( (x-xG)**2 +(y-yG)**2 +(z-zG)**2 ) + if minLength == -1 : + minLenght = length + minGroup = foundMultipleGroups[member] + minGroupIndex = member + continue + if minLength > length: + + minlegth = length + minGroup = foundMultipleGroups[member] + minGroupIndex = member +# +# the mighty try catch block ! +# +# + try: + xG = firstTubeDx.getGroupes()[minGroup][0][0] + oneMember = 0 + except TypeError: + xG = foundMultipleGroupsXYZ[minGroupIndex][0] + yG = foundMultipleGroupsXYZ[minGroupIndex][1] + zG = foundMultipleGroupsXYZ[minGroupIndex][2] + # newGroupes[minGroup].append([xG,yG,zG]) + newGroupes[minGroup].append([x,y,z]) + oneMember = 1 + newVoxelGroupGrid[x][y][z] = minGroup + if (oneMember == 0): +#### test this !!! +# print foundMultipleGroupsXYZ + #newGroupes[minGroup].append(foundMultipleGroupsXYZ[minGroupIndex]) + newGroupes[minGroup].append([x,y,z]) + newVoxelGroupGrid[x][y][z] = minGroup +#### this passage could be buggy either extend or append .... + + + + + + elif len(foundMultipleGroups) == 1: #### only one group was found + + try: + xG = firstTubeDx.getGroupes()[foundMultipleGroupsXYZ[0]][0][0] + oneMember = 0 + except TypeError: +# print foundMultipleGroupsXYZ + xG = foundMultipleGroupsXYZ[0][0] + yG = foundMultipleGroupsXYZ[0][1] + zG = foundMultipleGroupsXYZ[0][2] +# newGroupes[foundMultipleGroups[0]].append([xG,yG,zG]) + newGroupes[foundMultipleGroups[0]].append([x,y,z]) + oneMember = 1 + newVoxelGroupGrid[x][y][z] = foundMultipleGroupsXYZ[0] + if (oneMember == 0): +#### test this !!! + #newGroupes[minGroup].append(foundMultipleGroupsXYZ) + newGroupes[minGroup].append([x,y,z]) + newVoxelGroupGrid[x][y][z] = minGroup +#### this passage could be buggy either extend or append .... + outputLine = int ((x*100)/(len (self.__cavitysDxObject.getDensity()))) + print str(outputLine)+" %\r", + sys.stdout.flush() + len (self.__cavitysDxObject.getDensity()) + + print "Done" + newCavitiesTubeDx = TubeDx ( + openDxObject =self.__cavitysDxObject , # real densities + filename = self.__tubeDx.getFilename(), + proteinObject= self.__tubeDx.getProteinObject(), + waterObject= self.__tubeDx.getWaterObject(), + scanMethod = self.__tubeDx.getScanMethod(), + scanned=self.__tubeDx.getScanned(), + grouped='Yes', + groupes= newGroupes, + VoxelGroupGrid = newVoxelGroupGrid, + minDiameter = self.__tubeDx.getMinDiameter(), + protThreshold = self.__tubeDx.getProtThreshold(), + solventThreshold =self.__tubeDx.getSolventThreshold(), + version = self.__version, + protFile =self.__tubeDx.getProtFile(), + solvFile = self.__tubeDx.getSolvFile(), + filterApplied = self.__tubeDx.getFilterApplied(), + ) + + + return newCavitiesTubeDx + + + def joinByDistanceMajorGroup(self, firstTubeDx, secondOpenDx, distCutoff = 8): + +# first group second group ..... +#[ [x][y][z] [x][y][z] [x][y][z] ] + newGroupes = firstTubeDx.getGroupes() + +#[x][y][z] = groupID + newVoxelGroupGrid = firstTubeDx.getVoxelGroupGrid() + +#[x][y][z] = distance to a major group before + newVoxelGroupGridMinDist = [[[-1 for zR in range(len (self.__cavitysDxObject.getDensity()[0][0]))] for yR in range(len (self.__cavitysDxObject.getDensity()[0]))] for xR in range(len (self.__cavitysDxObject.getDensity()))] + i = 0 + groupLength = len (firstTubeDx.getGroupes()) + print "Cavities to process",groupLength + for groupID in firstTubeDx.getGroupes(): + + try: + xG = groupID[0][0] + oneMember = 0 + except TypeError: + oneMember = 1 + + + if oneMember: + size = 1 + while (size <= distCutoff): + negLen = size -size -size + posLen = size + + for xCurrent in range (negLen, (posLen+1)): + for yCurrent in range (negLen, (posLen+1)): + for zCurrent in range (negLen, (posLen+1)): + if xCurrent == 0 and yCurrent == 0 and zCurrent == 0: + continue + if firstTubeDx.getDxObject().getDensity()[xCurrent][yCurrent][zCurrent] > 0: # voxel is member of the "major" (first) selection + continue + xG = groupID[0] + yG = groupID[1] + zG = groupID[2] + + +# try : +# a = newVoxelGroupGrid[xCurrent+xG][yCurrent+yG][zCurrent+zG] +# if a == i: +# continue +# except IndexError: +# continue + + if (xG+xCurrent < len(firstTubeDx.getDxObject().getDensity()) and + yG+yCurrent < len (firstTubeDx.getDxObject().getDensity()[0]) and + zG+zCurrent < len( firstTubeDx.getDxObject().getDensity()[0][0]) and # checking if voxel is in array ... ... index out of bound error ... + xG+xCurrent >= 0 and yG+yCurrent >= 0 and zG+zCurrent >= 0 and + self.__tubeDx.getDxObject().getDensity()[xG+xCurrent][yG+yCurrent][zG+zCurrent] >= self.__tubeDx.getSolventThreshold() + ): +# x = xCurrent + xG +# xCurrent & yCurrent & zCurrent contain already delta x,y,z .... + length = math.sqrt( (xCurrent)**2 +(yCurrent)**2 +(zCurrent)**2 ) + + if ((length < newVoxelGroupGridMinDist[xCurrent+xG][yCurrent+yG][zCurrent+zG] and newVoxelGroupGridMinDist[xCurrent+xG][yCurrent+yG][zCurrent+zG] != -1) + or newVoxelGroupGridMinDist[xCurrent+xG][yCurrent+yG][zCurrent+zG] == -1): + newVoxelGroupGridMinDist[xCurrent+xG][yCurrent+yG][zCurrent+zG] = length + newVoxelGroupGrid[xCurrent+xG][yCurrent+yG][zCurrent+zG] = i + + + if oneMember == 0: + memberCount = 0 + currentGroupLength = len (groupID) + for member in groupID: + + size = 1 + while (size <= distCutoff): + negLen = size -size -size + posLen = size + + for xCurrent in range (negLen, (posLen+1)): + for yCurrent in range (negLen, (posLen+1)): + for zCurrent in range (negLen, (posLen+1)): + if xCurrent == 0 and yCurrent == 0 and zCurrent == 0: + continue +# if firstTubeDx.getDxObject().getDensity()[xCurrent][yCurrent][zCurrent] > 0: # voxel is member of the "major" (first) selection +# continue + + # only voxels above previosly defined threshold !! + xG = member[0] + yG = member[1] + zG = member[2] + +# try : +# a = newVoxelGroupGrid[xCurrent+xG][yCurrent+yG][zCurrent+zG] +# if a == i: +# continue +# except IndexError:# + #continue + + if (xG+xCurrent < len (firstTubeDx.getDxObject().getDensity()) and + yG+yCurrent < len (firstTubeDx.getDxObject().getDensity()[0]) and + zG+zCurrent < len (firstTubeDx.getDxObject().getDensity()[0][0]) and # checking if voxel is in array ... ... index out of bound error ... + xG+xCurrent >= 0 and yG+yCurrent >= 0 and zG+zCurrent >= 0 and + self.__tubeDx.getDxObject().getDensity()[xG+xCurrent][yG+yCurrent][zG+zCurrent] >= self.__tubeDx.getSolventThreshold() + ): + + # x = xCurrent + xG + # xCurrent & yCurrent & zCurrent contain already delta x,y,z .... + length = math.sqrt( (xCurrent)**2 +(yCurrent)**2 +(zCurrent)**2 ) +# print "yo" + if ((length < newVoxelGroupGridMinDist[xCurrent+xG][yCurrent+yG][zCurrent+zG] and newVoxelGroupGridMinDist[xCurrent+xG][yCurrent+yG][zCurrent+zG] != -1) + or newVoxelGroupGridMinDist[xCurrent+xG][yCurrent+yG][zCurrent+zG] == -1): + newVoxelGroupGridMinDist[xCurrent+xG][yCurrent+yG][zCurrent+zG] = length + newVoxelGroupGrid[xCurrent+xG][yCurrent+yG][zCurrent+zG] = i + size += 1 + + memberCount += 1 +# outputLine = int ((memberCount*100)/(currentGroupLength)) +# print str(outputLine)+" %\r", + print "Analysing member",memberCount,"/",currentGroupLength,"\r", + sys.stdout.flush() + + + +# outputLine = int ((i*100)/(len (firstTubeDx.getGroupes()))) +# print str(outputLine)+" %\r", + print "\nProcessed cavity",i + i += 1 + + + + for x in range (len (newVoxelGroupGridMinDist)): + for y in range (len (newVoxelGroupGridMinDist[0])): + for z in range (len (newVoxelGroupGridMinDist[0][0])): + if newVoxelGroupGridMinDist[x][y][z] >= self.__tubeDx.getSolventThreshold() : + newGroupes[newVoxelGroupGrid[x][y][z]].append([x,y,z]) + + + newCavitiesTubeDx = TubeDx ( + openDxObject =self.__cavitysDxObject , # real densities + filename = self.__tubeDx.getFilename(), + proteinObject= self.__tubeDx.getProteinObject(), + waterObject= self.__tubeDx.getWaterObject(), + scanMethod = self.__tubeDx.getScanMethod(), + scanned=self.__tubeDx.getScanned(), + grouped='Yes', + groupes= newGroupes, + VoxelGroupGrid = newVoxelGroupGrid, + minDiameter = self.__tubeDx.getMinDiameter(), + protThreshold = self.__tubeDx.getProtThreshold(), + solventThreshold =self.__tubeDx.getSolventThreshold(), + version = self.__version, + protFile =self.__tubeDx.getProtFile(), + solvFile = self.__tubeDx.getSolvFile(), + filterApplied = self.__tubeDx.getFilterApplied(), + ) + + + return newCavitiesTubeDx + +### +# This routine regroups cavities, each previously defined will be regrouped by grouptubes.py +# Regrouping is an advantage for join by distance routines cause theire results contain where not all voxels are connected +### + + + def reGroup(self, tubeDx): + print "Regrouping" +#groupes +#[[x][y][z][x][y][z]] + groupes =[] +#voxelGroupGrid +#[x][y][z] = group + densities = tubeDx.getDxObject().getDensity() + voxelGroupGrid = [[[-1 for zR in range(len (densities[0][0]))] for yR in range(len (densities[0]))] for xR in range(len (densities))] + for i in range (len(tubeDx.getGroupes())): + print "Group",i,"/",len(tubeDx.getGroupes()) + + currentGroupDensities = [[[0 for zR in range(len (densities[0][0]))] for yR in range(len (densities[0]))] for xR in range(len (densities))] + try: + xG = tubeDx.getGroupes()[i][0][0] + oneMember = 0 + except TypeError: + oneMember = 1 + if oneMember: + groupes.append(tubeDx.getGroupes()[i]) + voxelGroupGrid[tubeDx.getGroupes()[i][0]][tubeDx.getGroupes()[i][1]][tubeDx.getGroupes()[i][2]] = len(groupes) + if not oneMember: + oldLen = len (groupes) + + for member in tubeDx.getGroupes()[i]: +# x = member[0] +# y = member[1] +# z = member[2] + currentGroupDensities[member[0]][member[1]][member[2]] = tubeDx.getDxObject().getDensity()[member[0]][member[1]][member[2]] + + currentOpenDxObject = BuildOpenDx (density = currentGroupDensities, + stepSize = tubeDx.getDxObject().getStepsize() , + origin = tubeDx.getDxObject().getOrigin(), + dimention = tubeDx.getDxObject().getDimention() + ) + groupThis = GroupTubes (currentOpenDxObject) + groupThis.findGroups() + groupes.extend(groupThis.getTubeGroup()) + newLen = len (groupes) + if newLen - oldLen > 1: +# print "huh" + for groupIndex in range (oldLen, newLen-1): + for member in groupes[groupIndex]: +# print member,member[0],member[1],member[2] + voxelGroupGrid[member[0]][member[1]][member[2]] = groupIndex + else: + for member in groupes[newLen-1]: +# print member,member[0] + voxelGroupGrid[member[0]][member[1]][member[2]] = newLen-1 + + print len (groupes) +### grouping the rest ... :) + print "Grouping none core cavities" +# dxObject = tubeDx.getDxObject() + dxObject = self.__cavitysDxObject + currentGroupDensities = [[[0 for zR in range(len (densities[0][0]))] for yR in range(len (densities[0]))] for xR in range(len (densities))] + for x in range (len (dxObject.getDensity())): + for y in range (len (dxObject.getDensity()[0])): + for z in range (len (dxObject.getDensity()[0][0])): + if voxelGroupGrid[x][y][z] == -1 and dxObject.getDensity()[x][y][z] > 0: +# print "hu" + currentGroupDensities[x][y][z] = dxObject.getDensity()[x][y][z] + + currentOpenDxObject = BuildOpenDx (density = currentGroupDensities, + stepSize = tubeDx.getDxObject().getStepsize() , + origin = tubeDx.getDxObject().getOrigin(), + dimention = tubeDx.getDxObject().getDimention() + ) + oldLen = len (groupes) + groupThis = GroupTubes (currentOpenDxObject) + groupThis.findGroups() + groupes.extend(groupThis.getTubeGroup()) + newLen = len (groupes) + if newLen - oldLen > 1: +# print "huh" + for groupIndex in range (oldLen, newLen-1): + for member in groupes[groupIndex]: +# print member,member[0],member[1],member[2] + voxelGroupGrid[member[0]][member[1]][member[2]] = groupIndex + else: + for member in groupes[newLen-1]: +# print member,member[0] + voxelGroupGrid[member[0]][member[1]][member[2]] = newLen-1 + + + + print len (groupes) +# print groupes + newCavitiesTubeDx = TubeDx ( + openDxObject = tubeDx.getDxObject() , # real densities + filename = tubeDx.getFilename(), + proteinObject= tubeDx.getProteinObject(), + waterObject= tubeDx.getWaterObject(), + scanMethod = tubeDx.getScanMethod(), + scanned=tubeDx.getScanned(), + grouped='Yes', + groupes= groupes, + VoxelGroupGrid = voxelGroupGrid, + minDiameter = tubeDx.getMinDiameter(), + protThreshold = tubeDx.getProtThreshold(), + solventThreshold =tubeDx.getSolventThreshold(), + version = self.__version, + protFile =tubeDx.getProtFile(), + solvFile = tubeDx.getSolvFile(), + filterApplied = tubeDx.getFilterApplied(), + ) + return newCavitiesTubeDx + + +#### +# +# this is the "main" method for neighborhood driven postprocessing .... +# + + + def postGroup_Neighbor(self, threshold, cavityID, edgeLength = 2 , distCutoff = 8): + + neighborOpenDx = BuildOpenDx (density = self.calcNeighbors (cavityID = cavityID, edgeLength = edgeLength), + stepSize = self.__cavitysDxObject.getStepsize() , + origin = self.__cavitysDxObject.getOrigin(), + dimention = self.__cavitysDxObject.getDimention() + ) + print "Deviding cavity",cavityID,"in subcavites. Threshold =",threshold + densityArray = self.separateByValue (threshold, neighborOpenDx) + + firstOpenDx = BuildOpenDx (density = densityArray[0], + stepSize = self.__cavitysDxObject.getStepsize() , + origin = self.__cavitysDxObject.getOrigin(), + dimention = self.__cavitysDxObject.getDimention() + ) + secondOpenDx = BuildOpenDx (density = densityArray[1], + stepSize = self.__cavitysDxObject.getStepsize() , + origin = self.__cavitysDxObject.getOrigin(), + dimention = self.__cavitysDxObject.getDimention() + ) +# above threshold ! + + + firstGrouped = GroupTubes(firstOpenDx) +# print densityArray[0] + firstGrouped.findGroups() + + + tmpTubeDx = TubeDx ( + openDxObject =firstOpenDx , # <<<<----- neighbor values + filename = self.__tubeDx.getFilename(), + proteinObject= self.__tubeDx.getProteinObject(), + waterObject= self.__tubeDx.getWaterObject(), + scanMethod = self.__tubeDx.getScanMethod(), + scanned=self.__tubeDx.getScanned(), + grouped='Yes', + groupes= firstGrouped.getTubeGroup(), + VoxelGroupGrid = firstGrouped.getVoxelGroupGrid(), + minDiameter = self.__tubeDx.getMinDiameter(), + protThreshold = self.__tubeDx.getProtThreshold(), + solventThreshold =self.__tubeDx.getSolventThreshold(), + version = self.__version, + protFile =self.__tubeDx.getProtFile(), + solvFile = self.__tubeDx.getSolvFile(), + filterApplied = self.__tubeDx.getFilterApplied(), + ) + +# +# test ... lemme see if this can be written into pdb +# + print "### Debugging ###" + test = WritePDB (tmpTubeDx) + test.writeGroups(filename = 'core_cavities_after_filter.pdb', groupArray = firstGrouped.getTubeGroup()) + print "###/Debugging ###" + + + + + firstGroupedTubeDx = TubeDx ( + openDxObject =firstOpenDx , # <<<<----- neighbor values + filename = self.__tubeDx.getFilename(), + proteinObject= self.__tubeDx.getProteinObject(), + waterObject= self.__tubeDx.getWaterObject(), + scanMethod = self.__tubeDx.getScanMethod(), + scanned=self.__tubeDx.getScanned(), + grouped='Yes', + groupes= firstGrouped.getTubeGroup(), + VoxelGroupGrid = firstGrouped.getGroups(), + minDiameter = self.__tubeDx.getMinDiameter(), + protThreshold = self.__tubeDx.getProtThreshold(), + solventThreshold =self.__tubeDx.getSolventThreshold(), + version = self.__version, + protFile =self.__tubeDx.getProtFile(), + solvFile = self.__tubeDx.getSolvFile(), + filterApplied = self.__tubeDx.getFilterApplied(), + ) + print "Separating ... " +# return self.joinByDistance (firstTubeDx = firstGroupedTubeDx, secondOpenDx = secondOpenDx , distCutoff = 3) +### be carefull hard coded distance cutoff !!! +## +# this makes sence for acrb ... in future the distCutoff should be a cmd parameter ... +# + joinedTubeDx = self.joinByDistanceMajorGroup (firstTubeDx = firstGroupedTubeDx, secondOpenDx = secondOpenDx , distCutoff = distCutoff) + + return self.reGroup (tubeDx = joinedTubeDx) + + diff --git a/lib/postgroup.pyc b/lib/postgroup.pyc new file mode 100644 index 0000000..ae4b0f2 Binary files /dev/null and b/lib/postgroup.pyc differ diff --git a/lib/readopendx.py b/lib/readopendx.py new file mode 100644 index 0000000..b433311 --- /dev/null +++ b/lib/readopendx.py @@ -0,0 +1,163 @@ +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +###### +# v. 0.083 +#- floats as stepsize are now supported +# +# v. 0.066 +#- setting density to 0 if a negative density is read +#- removed bug if last row contains fever than 3 entries +###### +import re +from buildopendx import * +class ReadOpenDx(): + __dxfile =() + __origin =[] + __transformMatrix = [] + __size = [] + __stepSize = 0 + def __init__ (self, dxFile): + if (not dxFile): + print "Error no intputfile" +# exit + self.__dxFile = dxFile + + +# For details OpenDx fileformat visit opendx.org +# and http://apbs.wustl.edu/MediaWiki/index.php/OpenDX_scalar_data_format +# http://www.poissonboltzmann.org/file-formats/mesh-and-data-formats/opendx-scalar-data + + def importDx(self=None): + dxFile = self.__dxFile + dxFH = open(dxFile,"r") + dxlines = dxFH.readlines() + transformMatrix = [] + x=0 + y=0 + z=0 + row = 1 + for line in dxlines: + regExp = re.compile ('object 1') + if (regExp.match(line)): + self.__size = re.findall('\d+', line) + del self.__size[0] + for i in range(3): + self.__size[i] = int(self.__size[i]) + continue + + if (re.match('origin',line)): + self.__origin = re.findall('(-?[0-9]*\.?-?[0-9]+)', line) + self.__origin[0] = float(self.__origin[0]) + self.__origin[1] = float(self.__origin[1]) + self.__origin[2] = float(self.__origin[2]) + print 'Origin \nx: '+str(self.__origin[0])+'\ny: '+str(self.__origin[1])+'\nz: '+str(self.__origin[2]) + continue + +# transformMatrix: +# transformMatrix[0] transformation to get to middle of the continue boxes x-value +# transformMatrix[1] transformation to get to middle of the continue boxes y-value etc ... + if (re.match('delta',line)): + if self.__stepSize > 0: + continue +# self.__transformMatrix.append (re.findall('(-?[0-9]*\.?-?[0-9]+)', line)) +# self.__transformMatrix.append (re.findall('\d*\.*\d+', line)) + self.__stepSize = (re.findall('\d*\.*\d+', line)) [0] + print "Stepsize: " + str(self.__stepSize) + continue + if (re.match('^-*[0-9]', line)): +# initialysing density array [x][y][z] = density + if (x == 0 and y == 0 and z == 0): + density = [[[0 for zR in range(self.__size[2])] for yR in range(self.__size[1])] for xR in range(self.__size[0])] +# 3 densitys in one row ... z is increasing fastest y medium x slowest +# for openDx file format read : http://apbs.wustl.edu/MediaWiki/index.php/OpenDX_scalar_data_format +# Remember x y and z are multipliers for transformmatrix to get to the middle of the continue box + densityRow = (re.split('\s', line)) + densityRow.pop() + for i in range (len(densityRow)): + if (densityRow[i] == '' ): + densityRow.pop() + continue + if float(densityRow[i]) < 0: + densityRow[i] = 0 + if (z >= self.__size[2]): + z = 0 + y += 1 + if (y >= self.__size[1]): + y = 0 + x += 1 + if (x >= self.__size[0]): + x = 0 + density[x][y][z] = float(densityRow[0]) + z += 1 + + if (x >= self.__size[0] -1 and y >= self.__size[1] - 1 and z >= self.__size[2] -1): + break + if (z >= self.__size[2]): + z = 0 + y += 1 + if (y >= self.__size[1]): + y = 0 + x += 1 + if (x >= self.__size[0]): + x = 0 + density[x][y][z] = float(densityRow[1]) + z += 1 + + if (x >= self.__size[0] -1 and y >= self.__size[1] - 1 and z >= self.__size[2] -1): + break + if (z >= self.__size[2]): + z = 0 + y += 1 + if (y >= self.__size[1]): + y = 0 + x += 1 + if (x >= self.__size[0]): + x = 0 + density[x][y][z] = float(densityRow[2]) + z += 1 + outputLine = int ((row*100)/(( (self.__size[0] *self.__size[1]*self.__size[2])/3) )) + print str(outputLine)+" %\b\b\b\b\b", + row += 1 + dxFH.close() + self.__density = density + print "Done" + def getOrigin(self = None): + return self.__origin + def getSize(self=None): + tmp = int(int(self.__size[0])* int(self.__size[1])*int(self.__size[2])) + return tmp + def getDimention(self=None): + return self.__size + def getStepsize(self=None): + return float(self.__stepSize) + def getDensity(self=None): + return self.__density + + def setDxfile(self,dxFile): + if (not dxFile): + print "Error no intputfile" + exit + self.__dxFile = dxFile + self.__converted = () + self.__origin =() + self.__transformMatrix = [] + self.__size = () + def getFilename(self): + return self.__dxFile + def getDxObject(self=None): + object = BuildOpenDx (self.__density, float(self.__stepSize), self.__origin,self.__size) + return object diff --git a/lib/readopendx.pyc b/lib/readopendx.pyc new file mode 100644 index 0000000..8b3d747 Binary files /dev/null and b/lib/readopendx.pyc differ diff --git a/lib/readpdb.py b/lib/readpdb.py new file mode 100644 index 0000000..ff24189 --- /dev/null +++ b/lib/readpdb.py @@ -0,0 +1,232 @@ +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +### +# +#v0.27 +#+Added error message for unsupported pdb files. Example header is now printed to cli +# +#v0.12 +# -coordinates are extrated now absolute postions +# + + + +import re +import sys +from tubedx import * +from buildopendx import * + +class ReadPDB (): + __scanned = 'No' + __filename="None" + __scanType = 'None' + __minSolvDens ='' + __minProtDens ='' + __minDiameter = '' + __protFile = '' + __solvFile = '' + __version = 0 + __filter = '' + __grouped = 'No' + def __init__ (self, pdbFile): + self.__pdbFile = pdbFile + def readPDB (self): + pdbFH = open(self.__pdbFile,"r") + origin = ['x','y','z'] + maxCoord = ['x','y','z'] + dimention = [0,0,0] + groupArray = [] + lastGroup = 0 + for line in pdbFH.readlines(): +#getting dxTuber settings + if (re.match('^dxTuber\s+Version',line)): + self.__version = line.rsplit()[2] + print 'Created with Version ' + self.__version + if (re.match('^dxTuber\s+OpenDXFile',line)): + self.__dxFile = line.rsplit()[2] +# print line.rsplit()[2] + + if (re.match('^dxTuber\s+ProteinFile',line)): + self.__protFile = line.rsplit()[2] + print 'Protein file ' + self.__protFile + if (re.match('^dxTuber\s+SolventFile',line)): + self.__solvFile = line.rsplit()[2] + print 'Solvent file ' + self.__solvFile + if (re.match('^dxTuber\s+StepSize',line)): + self.__stepSize = line.rsplit()[2] + print 'Stepsize ' + self.__stepSize + if (re.match('^dxTuber\s+ScanType',line)): + self.__scanType = line.rsplit()[2] + print 'Scan type ' + self.__scanType + if (re.match('^dxTuber\s+MinDiameter',line)): + self.__minDiameter = line.rsplit()[2] + print 'Min diameter ' + self.__minDiameter + if (re.match('^dxTuber\s+MinProtDens',line)): + self.__minProtDens = line.rsplit()[2] + print 'Min protein density ' +self.__minProtDens + if (re.match('^dxTuber\s+MinSolvDens',line)): + self.__minSolvDens= line.rsplit()[2] + print 'Min solvent density ' + self.__minSolvDens + if (re.match('^dxTuber\s+Filter',line)): +#filters are seperated by ', ' + tempLine = line.rstrip('\n') + self.__filter= tempLine.split(None,2)[2] + print 'Filter applied ' + self.__filter +#getting the origin +# x_coord = line[31:38] +# y_coord= line[39:46] +# z_coord= line[47:54] + + + if (re.match('^ATOM',line)): +# print 'x: ' + line[31:38] +# print 'y: ' + line[39:46] +# print 'z: ' + line[47:54] + if origin[0] == 'x': + origin[0] = float (line[31:38]) + origin[1] = float (line[39:46]) + origin[2] = float (line[47:54]) + + if (origin[1] > float (line[39:46])): + origin[1] = float (line[39:46]) + if (origin[2] > float (line[47:54])): + origin[2] = float (line[47:54]) +#getting size of the grid + if (maxCoord[0] == 'x'): + maxCoord[0] = float (line[31:38]) + maxCoord[1] = float (line[39:46]) + maxCoord[2] = float (line[47:54]) + + if (maxCoord[0] < float (line[31:38])): + maxCoord[0] = float (line[31:38]) + if (maxCoord[1] < float (line[39:46])): + maxCoord[1] = float (line[39:46]) + if (maxCoord[2] < float (line[47:54])): + maxCoord[2] = float (line[47:54]) + pdbFH.close() + if self.__version == 0: + print 'PDB file not supported \nOnly previously processed files can be read.' +# print 'Please add an "dxTuber" header to your file (e.g.): \ndxTuber Version 0.25\ndxTuber OpenDXFile None\ndxTuber ProteinFile 110_protein.dx\ndxTuber SolventFile 110_water.dx\ndxTuber StepSize 1.0\ndxTuber ScanType 2D\ndxTuber MinDiameter 0\ndxTuber MinProtDens 0.005\ndxTuber MinSolvDens 0.0\ndxTuber Filter None\ndxTuber Grouped No\n' + print """ +Please add an "dxTuber" header to your file (e.g.): +dxTuber Version 0.25 +dxTuber OpenDXFile None +dxTuber ProteinFile 110_protein.dx +dxTuber SolventFile 110_water.dx +dxTuber StepSize 1.0 +dxTuber ScanType 2D +dxTuber MinDiameter 0 +dxTuber MinProtDens 0.005 +dxTuber MinSolvDens 0.0 +dxTuber Filter None +dxTuber Grouped No""" +# + sys.exit() +# return (1) + +# grid should not end directly after a cavity to make grouping possible + maxCoord[0] = maxCoord[0] +2 + maxCoord[1] = maxCoord[1] +2 + maxCoord[2] = maxCoord[2] +2 + + origin[0] = origin[0] -1 + origin[1] = origin[1] -1 + origin[2] = origin[2] -1 + + dimention[0] = int(maxCoord[0] - origin[0]) + dimention[1] = int(maxCoord[1] - origin[1]) + dimention[2] = int(maxCoord[2] - origin[2]) + print 'Origin:' + print origin + print 'Dimension:' + print dimention + density = [[[0 for zR in range(dimention[2])] for yR in range(dimention[1])] for xR in range(dimention[0])] +#groupGrid[x][y][z] = group array + groupGrid = [[[-1 for zR in range(dimention[2])] for yR in range(dimention[1])] for xR in range(dimention[0])] +# fetching densities + pdbFH = open(self.__pdbFile,"r") + + for line in pdbFH.readlines(): + if (re.match('^ATOM',line)): + x = int(float(line[31:38]) - origin[0]) + y = int(float(line[39:46]) - origin[1]) + z = int(float(line[47:54]) - origin[2]) + density[x][y][z] = float (line[61:66]) + if re.search('[0-9]+', line[13:16]): + self.__grouped = 'Yes' + if groupArray == []: + tmp =[x,y,z] + groupArray.append([tmp]) + elif int(line[13:16]) == lastGroup: + groupArray[lastGroup].append ([x,y,z]) + self.__grouped = 'Yes' + if int(line[13:16]) != lastGroup: + tmp = [([x,y,z])] + groupArray.append (tmp) + lastGroup += 1 +#groupGrid[x][y][z] = group array + if groupGrid[x][y][z] != -1: + groupGrid[x][y][z] +=' '+ str(int(line[13:16])) + else: + groupGrid[x][y][z] = str(int(line[13:16])) + pdbFH.close() + self.__density = density + self.__origin = origin + self.__size = dimention + if self.__scanType != 'None': + self.scanned = 'Yes' + try: + self.__dxFile + except: + self.__dxFile = 'None' + if self.__grouped == 'Yes': + self.__tubeDxObject = TubeDx ( self.getDxObject(), + filename=self.__dxFile, + scanMethod = self.__scanType, + scanned = self.__scanType, + solventThreshold = self.__minSolvDens , + protThreshold = self.__minProtDens, + minDiameter = self.__minDiameter, + protFile = self.__protFile, + solvFile = self.__solvFile, + version = self.__version, + filterApplied = self.__filter, + groupes = groupArray, + VoxelGroupGrid = groupGrid, + grouped = 'Yes' + ) + else: + self.__tubeDxObject = TubeDx ( self.getDxObject(), + filename= self.__dxFile, + scanMethod = self.__scanType, + scanned = self.__scanType, + solventThreshold = self.__minSolvDens , + protThreshold = self.__minProtDens, + minDiameter = self.__minDiameter, + protFile = self.__protFile, + solvFile = self.__solvFile, + version = self.__version, + filterApplied = self.__filter + ) + print "Finished importing: "+ self.__pdbFile + return (0) + def getDxObject(self=None): + object = BuildOpenDx (self.__density, self.__stepSize, self.__origin,self.__size) + return object + def getTubeDxObject (self): + return self.__tubeDxObject + diff --git a/lib/readpdb.pyc b/lib/readpdb.pyc new file mode 100644 index 0000000..e5cb8f2 Binary files /dev/null and b/lib/readpdb.pyc differ diff --git a/lib/settings.py b/lib/settings.py new file mode 100644 index 0000000..bb54c65 --- /dev/null +++ b/lib/settings.py @@ -0,0 +1,43 @@ +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +### + + + +# dxTuber v0.2 +# +new routine to parse settings from $HOME/.dxTuber/dxTuber.conf +# + +import os, re + +class Settings(): + + + def getVersion(self): + dxTuberConfFile = os.environ['HOME']+"/.dxTuber/dxTuber.conf" + conf = open (dxTuberConfFile,'r') + confLines = conf.readlines() + + for line in confLines: + if re.match ('^#',line): + continue + if re.match('VERSION',line): + tmpArray = re.split('\s',line,1) + tmp = re.sub('\"','',tmpArray[1]) + tmp = re.sub('\n','',tmp) + version = re.sub('^\s','',tmp) + return version diff --git a/lib/settings.pyc b/lib/settings.pyc new file mode 100644 index 0000000..eb8bcb7 Binary files /dev/null and b/lib/settings.pyc differ diff --git a/lib/tubedx.py b/lib/tubedx.py new file mode 100644 index 0000000..decce53 --- /dev/null +++ b/lib/tubedx.py @@ -0,0 +1,191 @@ +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +### +# +# v0.15 +# +supporting +# +dxObjectArray +# +groupesArray +# +VoxelGroupGridArray + + +class TubeDx(): + def __init__(self, + openDxObject, + filename = 'None', + proteinObject=0, + waterObject=0, + scanMethod = '-1', + scanned='No', + grouped='No', + groupes=None, + VoxelGroupGrid = None, + minDiameter = '0' , + protThreshold = '0', + solventThreshold ='0', + version = 'unknown', + protFile = 'None', + solvFile = 'None', + filterApplied = 'None', + +# following arrays could be used for filtering: +# +# position 0 contain the Voxels fullfilling the filter criterion openDxObjects / groupes / voxelgroupgrids +# position 1 contain the Voxels not fullfilling the filter criterion openDxObjects / groupes / voxelgroupgrids +# +# both posiotions will be grouped separately and written into PDB file +# + + +#since v0.15 containing multiple dxObjects +# +# minimum diameter scanning +# dxObjectArray[0] densities bigger than minDiameter threshold +# dxObjectArray[1] densities smaler than minDiameter threshold +# + dxObjectArray = None, + groupesArray = None, + VoxelGroupGridArray = None + ): + + self.__openDxObject = openDxObject + self.__filename = filename + self.__proteinObject = proteinObject + self.__waterObject = waterObject + self.__scanned = scanned # if the cavity is scanned types are '1D[xyz] , 2D , 3D .... + self.__grouped = grouped + self.__version = version + self.__scanMethod = str(scanMethod) + self.__protThreshold = float(protThreshold) + self.__solventThreshold = float(solventThreshold) + self.__minDiameter = minDiameter + self.__protFile = protFile + self.__solvFile = solvFile + self.__filterApplied = filterApplied +# Be aware of these different arrays ... +#groupes +#[[x][y][z][x][y][z]] + self.__groupes = groupes +#voxelGroupGrid +#[x][y][z] = group + self.__VoxelGroupGrid = VoxelGroupGrid + +# since v0.15: +# +# [0] Voxels fullfilling the filter criterion +# [1] Voxels not fullfilling the filter criterion + self.__dxObjectArray = dxObjectArray + self.__groupesArray = groupesArray + self.__VoxelGroupGridArray = VoxelGroupGridArray + + + + def getFilename(self): + return self.__filename + def getProtFile(self): + return self.__protFile + def getSolvFile(self): + return self.__solvFile + def setFilterApplied(self,filter): + if (self.__filterApplied == 'None'): + self.__filterApplied = filter + else: self.__filterApplied += ', '+ filter + def getFilterApplied(self): + return self.__filterApplied + + def setFilename(self,filename): + self.__filename = filename + def setWater (self,condition = 0): + self.__waterObject = condition + def getWater (self): + return self.__waterObject + def getVersion(self): + return self.__version + def getProtThreshold(self): + return self.__protThreshold + def getSolventThreshold(self): + return self.__solventThreshold + def getMinDiameter(self): + return self.__minDiameter + def setProtein (self,condition = 0): + self.__proteinObject = condition + def getProtein (self): + return self.__proteinObject + + def setScanned (self,condition = 0): + self.__scanned = condition + def setGrouped (self,condition = 0): + self.__grouped = condition + def setScanMethod (self,condition = 0): + self.__scanMethod = condition + def getScanMethod (self): + return self.__scanMethod +#groupes +#[[x][y][z][x][y][z]] + def setGroupes (self,condition = None): + self.__groupes = condition + def getGroupes (self): + return self.__groupes +#voxelGroupGrid +#[x][y][z] = group + def setVoxelGroupGrid (self,condition = None): + self.__VoxelGroupGrid = condition + def getVoxelGroupGrid (self): + return self.__VoxelGroupGrid + + + def getScanned (self): + return self.__scanned + def getGrouped (self): + return self.__grouped + + def getDxObject(self): + return self.__openDxObject + def setDxObject(self,dxObject): + self.__openDxObject = dxObject + + def setProteinObject(self,openDxObject): + self.__proteinObject = openDxObject + def getProteinObject(self): + return self.__proteinObject + def setWaterObject(self,openDxObject): + self.__waterObject = openDxObject + def getWaterObject(self): + return self.__waterObject + +##### new in v0.15 +# +# + + def getDxObjectArray(self): + return self.__dxObjectArray + +#array of array of groupes +#[ [[x][y][z][x][y][z]] ] <--- 0 and more like them 1 ----> [ [[x][y][z][x][y][z]] ] + def getGroupesArray(self): + return self.__groupesArray + def getVoxelGroupGridArray(self): + return self.__VoxelGroupGridArray + + def setDxObjectArray(self, condition): + self.__dxObjectArray = condition + def setGroupesArray(self,condition): + self.__groupesArray = condition + def setVoxelGroupGridArray(self,condition): + self.__VoxelGroupGridArray = condition + + diff --git a/lib/tubedx.pyc b/lib/tubedx.pyc new file mode 100644 index 0000000..4caddcb Binary files /dev/null and b/lib/tubedx.pyc differ diff --git a/lib/writedx.py b/lib/writedx.py new file mode 100644 index 0000000..d271afa --- /dev/null +++ b/lib/writedx.py @@ -0,0 +1,83 @@ +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. + +class WriteDx(): + __dxObject = () + def __init__(self, dxObject): + self.__dxObject = dxObject + def write(self,filename): + dxFH = open(filename,"w") + density = self.__dxObject.getDensity() + origin = self.__dxObject.getOrigin() + atomCount = len (density)*len (density[0])*len (density[0][0]) + triple = () + dxFH.write ("# Created via dxTuber\n") + dxFH.write ("object 1 class gridpositions counts %i %i %i\n" %(len (density),len (density[0]),len (density[0][0]))) + dxFH.write ("origin %.1f %.1f %.1f\n" %(origin[0],origin[1],origin[2])) +#writing transformation matrx + dxFH.write ("delta %s 0 0\n" %(self.__dxObject.getStepsize()) ) + dxFH.write ("delta 0 %s 0\n" %(self.__dxObject.getStepsize()) ) + dxFH.write ("delta 0 0 %s\n" %(self.__dxObject.getStepsize()) ) + dxFH.write ("object 2 class gridconnections counts %i %i %i\n" %(len (density),len (density[0]),len (density[0][0]))) + dxFH.write ("object 3 class array type double rank 0 items %i data follows\n" %(atomCount)) + +#header complete + triple= [] + for xCount in range (len (density)): + for yCount in range (len (density[0])): + for zCount in range (len (density[0][0])): + triple.append (density[xCount][yCount][zCount]) + if (len (triple) == 3): + for i in range (len (triple)): + if i < 2: + if (triple[i] == 0): + dxFH.write ("0 ") + else: + dxFH.write ("%e " %(triple[i])) + else: + if (triple[i] == 0): + dxFH.write ("0") + else: + dxFH.write ("%e" %(triple[i])) + dxFH.write ("\n") + triple = [] + +#if one entry missing + if (len (triple) == 1 ): + if (triple[0] == 0): + dxFH.write ("0") + else: + dxFH.write ("%e" %(triple[0])) + dxFH.write ("\n") + +#if 2 entries missing + if (len (triple) == 2 ): + if (triple[0] == 0): + dxFH.write ("0 ") + else: + dxFH.write ("%e " %(triple[0])) + if (triple[1] == 0): + dxFH.write ("0") + else: + dxFH.write ("%e" %(triple[1])) + dxFH.write ("\n") + + + dxFH.write ("\nobject \"density (protein) [A^-3]\" class field") + dxFH.close () + print "Saved OpenDX in "+str(filename) + diff --git a/lib/writedx.pyc b/lib/writedx.pyc new file mode 100644 index 0000000..3807d88 Binary files /dev/null and b/lib/writedx.pyc differ diff --git a/lib/writepdb.py b/lib/writepdb.py new file mode 100644 index 0000000..9b5ce54 --- /dev/null +++ b/lib/writepdb.py @@ -0,0 +1,252 @@ +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +# +### +# +# v0.2 +# -ungrouped cavites have group id '0' instead of 'DEN' -> analyzeCavity.py can process ungrouped cavities now +# +writeMultiGroupeArray (implemented for mdk (minimum distance keep)) +# -groupe id's are not printed anymore to cmd // stop spamming +# -densities > 99 are now supported // densities > 99 6.2f > 999 6.1f > 9999 6.0f > 999999 strange things could will happen :> +# +# v0.09 +# +# +Grouped Yes / No field in header +# +saved as +# +# v0.085 +# +'cryst' entry to enable pbc in molecular viewer +# +# v0.082 +# +dxTuber settings are saved in pdb files +# +# v0.08 +# +user defined solvent min densities are filtered already during scanningmethod ! +# +# +## +class WritePDB (): + __density=[] + __origin =[] + __stepSize = () + __origin = () + __dimention = [] + def __init__(self, tubeDx): + dxObject = tubeDx.getDxObject() + self.__density = dxObject.getDensity() + self.__stepSize = dxObject.getStepsize() + self.__origin = dxObject.getOrigin() + self.__tubeDx = tubeDx + def write (self, filename): + i=0 + pdbFH = open(filename,"w") + atomCount = 1 +#begin dxTuber parameter section + + pdbFH.write("dxTuber\tVersion\t\t"+self.__tubeDx.getVersion()+"\n") + pdbFH.write("dxTuber\tOpenDXFile\t"+self.__tubeDx.getFilename()+"\n") + pdbFH.write("dxTuber\tProteinFile\t"+self.__tubeDx.getProtFile()+"\n") + pdbFH.write("dxTuber\tSolventFile\t"+self.__tubeDx.getSolvFile()+"\n") + pdbFH.write("dxTuber\tStepSize\t"+str(self.__tubeDx.getDxObject().getStepsize())+"\n") + pdbFH.write("dxTuber\tScanType\t"+self.__tubeDx.getScanned()+"\n") + pdbFH.write("dxTuber\tMinDiameter\t"+self.__tubeDx.getMinDiameter()+"\n") + pdbFH.write("dxTuber\tMinProtDens\t"+str(self.__tubeDx.getProtThreshold())+"\n") + pdbFH.write("dxTuber\tMinSolvDens\t"+str(self.__tubeDx.getSolventThreshold())+"\n") + pdbFH.write("dxTuber\tFilter\t\t"+self.__tubeDx.getFilterApplied()+"\n") + pdbFH.write("dxTuber\tGrouped\t\tNo\n") +#CRYST1 98.959 98.084 186.817 90.00 90.00 90.00 P 1 1 + pdbFH.write ("CRYST1%9.3f%9.3f%9.3f 90.00 90.00 90.00 P 1 1\n" %(len (self.__density),len (self.__density[0]),len (self.__density[0][0]))) + pdbFH.write("\n") + + for xCount in range (len (self.__density)): + for yCount in range (len (self.__density[0])): + for zCount in range (len (self.__density[0][0])): +#For visual representation we considered only voxels exceeding a cutoff density +#of 0.015 H2O/A3, which is ~50% of the value for bulk water under standard conditions +#http://www.pubmedcentral.nih.gov/articlerender.fcgi?tool=pubmed&pubmedid=14747309 +# +#Dynamics of Water Molecules in the Bacteriorhodopsin Trimer in Explicit Lipid/Water Environment +#Christian Kandt, Juergen Schlitter, and Klaus Gerwert +#so here cutoff 0.01 H_2O was choosen +# not nessaccary here, user can define this via gui / cmd and is done by planesearch !! + if (self.__density[xCount][yCount][zCount] > 0): + x = float (xCount) * float(self.__stepSize) + float(self.__origin[0]) + y = float (yCount) * float(self.__stepSize) + float(self.__origin[1]) + z = float (zCount) * float(self.__stepSize) + float(self.__origin[2]) + if (atomCount > 99999): + atomCount = 1 +# pdbFH.write ("ATOM %5i DEN DEN H 1 %8.3f%8.3f%8.3f 1.00%6.3f D\n" %( atomCount,x, y, z, self.__density[xCount][yCount][zCount]) ) + if self.__density[xCount][yCount][zCount] < 100: + pdbFH.write ("ATOM %5i 0 DEN H 1 %8.3f%8.3f%8.3f 1.00%6.3f D\n" %( atomCount,x, y, z, self.__density[xCount][yCount][zCount]) ) + elif self.__density[xCount][yCount][zCount] < 1000: + pdbFH.write ("ATOM %5i 0 DEN H 1 %8.3f%8.3f%8.3f 1.00%6.2f D\n" %( atomCount,x, y, z, self.__density[xCount][yCount][zCount]) ) + elif self.__density[xCount][yCount][zCount] < 10000: + pdbFH.write ("ATOM %5i 0 DEN H 1 %8.3f%8.3f%8.3f 1.00%6.1f D\n" %( atomCount,x, y, z, self.__density[xCount][yCount][zCount]) ) + else: + pdbFH.write ("ATOM %5i 0 DEN H 1 %8.3f%8.3f%8.3f 1.00%6.0f D\n" %( atomCount,x, y, z, self.__density[xCount][yCount][zCount]) ) + atomCount += 1 + pdbFH.write("TER\n") + pdbFH.write("END\n") + + pdbFH.close() + print "saved "+filename + +## group array +#groupes +#[[x][y][z][x][y][z]] + def writeGroups(self, groupArray, filename): + atomCount = 1 + chain =0 + maxChain = len (groupArray) + pdbFH = open(filename, "w") + + pdbFH.write("dxTuber\tVersion\t\t"+self.__tubeDx.getVersion()+"\n") + pdbFH.write("dxTuber\tOpenDXFile\t"+self.__tubeDx.getFilename()+"\n") + pdbFH.write("dxTuber\tProteinFile\t"+self.__tubeDx.getProtFile()+"\n") + pdbFH.write("dxTuber\tSolventFile\t"+self.__tubeDx.getSolvFile()+"\n") + pdbFH.write("dxTuber\tStepSize\t"+str(self.__tubeDx.getDxObject().getStepsize())+"\n") + pdbFH.write("dxTuber\tScanType\t"+self.__tubeDx.getScanned()+"\n") + pdbFH.write("dxTuber\tMinDiameter\t"+self.__tubeDx.getMinDiameter()+"\n") + pdbFH.write("dxTuber\tMinProtDens\t"+str(self.__tubeDx.getProtThreshold())+"\n") + pdbFH.write("dxTuber\tMinSolvDens\t"+str(self.__tubeDx.getSolventThreshold())+"\n") + pdbFH.write("dxTuber\tFilter\t\t"+self.__tubeDx.getFilterApplied()+"\n") + pdbFH.write("dxTuber\tGrouped\t\tYes\n") +#CRYST1 98.959 98.084 186.817 90.00 90.00 90.00 P 1 1 + pdbFH.write ("CRYST1%9.3f%9.3f%9.3f 90.00 90.00 90.00 P 1 1\n" %(len (self.__density),len (self.__density[0]),len (self.__density[0][0]))) + + pdbFH.write("\n") + + + while (chain < maxChain): +# print chain #// stop spamming ... ;) + for mCounter in range (len (groupArray[chain])): + if (atomCount > 99999): + atomCount = 1 + try: + x = float (float(groupArray[chain][mCounter][0]) * float(self.__stepSize)) + float(self.__origin[0]) + y = float (float(groupArray[chain][mCounter][1]) * float(self.__stepSize)) + float(self.__origin[1]) + z = float (float(groupArray[chain][mCounter][2]) * float(self.__stepSize)) + float(self.__origin[2]) +### scaling output format for densities ... +# v0.2 densities up to 999999 are now supported ... + if self.__density[groupArray[chain][mCounter][0]][groupArray[chain][mCounter][1]][groupArray[chain][mCounter][2]] < 100: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.3f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][mCounter][0]][groupArray[chain][mCounter][1]][groupArray[chain][mCounter][2]]) ) + elif self.__density[groupArray[chain][mCounter][0]][groupArray[chain][mCounter][1]][groupArray[chain][mCounter][2]] < 1000: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.2f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][mCounter][0]][groupArray[chain][mCounter][1]][groupArray[chain][mCounter][2]]) ) + elif self.__density[groupArray[chain][mCounter][0]][groupArray[chain][mCounter][1]][groupArray[chain][mCounter][2]] < 10000: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.1f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][mCounter][0]][groupArray[chain][mCounter][1]][groupArray[chain][mCounter][2]]) ) + else: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.0f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][mCounter][0]][groupArray[chain][mCounter][1]][groupArray[chain][mCounter][2]]) ) + + except TypeError: + x = float (float(groupArray[chain][0]) * float(self.__stepSize)) + float(self.__origin[0]) + y = float (float(groupArray[chain][1]) * float(self.__stepSize)) + float(self.__origin[1]) + z = float (float(groupArray[chain][2]) * float(self.__stepSize)) + float(self.__origin[2]) + + if self.__density[groupArray[chain][0]][groupArray[chain][1]][groupArray[chain][2]] < 100: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.3f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][0]][groupArray[chain][1]][groupArray[chain][2]]) ) + elif self.__density[groupArray[chain][0]][groupArray[chain][1]][groupArray[chain][2]] < 1000: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.2f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][0]][groupArray[chain][1]][groupArray[chain][2]]) ) + elif self.__density[groupArray[chain][0]][groupArray[chain][1]][groupArray[chain][2]] < 10000: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.1f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][0]][groupArray[chain][1]][groupArray[chain][2]]) ) + else: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.0f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][0]][groupArray[chain][1]][groupArray[chain][2]]) ) + + atomCount += 1 + chain += 1 + pdbFH.write("TER\n") + pdbFH.write("END\n") + + pdbFH.close + print "saved "+filename +## MultiGroupArray contains arrays of groups // [0] array one [1]array two .... +# +# Multiple groups may contain +# [0] groups of isv that fullfill a mindiameter threshold +# [1] groups of isv that do not fullfill a mindiameter threshold +# +## member member +# { [[x][y][z][x][y][z]] } { [[x][y][z][x][y][z]] } +# G R O U P Array +# + def writeMultiGroups(self, multiGroupArray, filename): + print "Saving multiple Groups" + atomCount = 1 + chain =0 #### GROAAAR this row took 1 week of bug fixing/searching m( + pdbFH = open(filename, "w") + + pdbFH.write("dxTuber\tVersion\t\t"+self.__tubeDx.getVersion()+"\n") + pdbFH.write("dxTuber\tOpenDXFile\t"+self.__tubeDx.getFilename()+"\n") + pdbFH.write("dxTuber\tProteinFile\t"+self.__tubeDx.getProtFile()+"\n") + pdbFH.write("dxTuber\tSolventFile\t"+self.__tubeDx.getSolvFile()+"\n") + pdbFH.write("dxTuber\tStepSize\t"+str(self.__tubeDx.getDxObject().getStepsize())+"\n") + pdbFH.write("dxTuber\tScanType\t"+self.__tubeDx.getScanned()+"\n") + pdbFH.write("dxTuber\tMinDiameter\t"+self.__tubeDx.getMinDiameter()+"\n") + pdbFH.write("dxTuber\tMinProtDens\t"+str(self.__tubeDx.getProtThreshold())+"\n") + pdbFH.write("dxTuber\tMinSolvDens\t"+str(self.__tubeDx.getSolventThreshold())+"\n") + pdbFH.write("dxTuber\tFilter\t\t"+self.__tubeDx.getFilterApplied()+"\n") + pdbFH.write("dxTuber\tGrouped\t\tYes\n") +#CRYST1 98.959 98.084 186.817 90.00 90.00 90.00 P 1 1 + pdbFH.write ("CRYST1%9.3f%9.3f%9.3f 90.00 90.00 90.00 P 1 1\n" %(len (self.__density),len (self.__density[0]),len (self.__density[0][0]))) + + pdbFH.write("\n") + + for groupArray in multiGroupArray: + chain =0 + print str(len(multiGroupArray))+" len multiGroupArray" + maxChain = len(groupArray) + while (chain < maxChain): +# print chain #// stop spamming ... ;) + for mCounter in range (len (groupArray[chain])): + if (atomCount > 99999): + atomCount = 1 + try: + x = float (float(groupArray[chain][mCounter][0]) * float(self.__stepSize)) + float(self.__origin[0]) + y = float (float(groupArray[chain][mCounter][1]) * float(self.__stepSize)) + float(self.__origin[1]) + z = float (float(groupArray[chain][mCounter][2]) * float(self.__stepSize)) + float(self.__origin[2]) + + if self.__density[groupArray[chain][mCounter][0]][groupArray[chain][mCounter][1]][groupArray[chain][mCounter][2]] < 100: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.3f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][mCounter][0]][groupArray[chain][mCounter][1]][groupArray[chain][mCounter][2]]) ) + elif self.__density[groupArray[chain][mCounter][0]][groupArray[chain][mCounter][1]][groupArray[chain][mCounter][2]] < 1000: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.2f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][mCounter][0]][groupArray[chain][mCounter][1]][groupArray[chain][mCounter][2]]) ) + elif self.__density[groupArray[chain][mCounter][0]][groupArray[chain][mCounter][1]][groupArray[chain][mCounter][2]] < 10000: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.1f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][mCounter][0]][groupArray[chain][mCounter][1]][groupArray[chain][mCounter][2]]) ) + else: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.0f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][mCounter][0]][groupArray[chain][mCounter][1]][groupArray[chain][mCounter][2]]) ) + + + except TypeError: + x = float (float(groupArray[chain][0]) * float(self.__stepSize)) + float(self.__origin[0]) + y = float (float(groupArray[chain][1]) * float(self.__stepSize)) + float(self.__origin[1]) + z = float (float(groupArray[chain][2]) * float(self.__stepSize)) + float(self.__origin[2]) + + if self.__density[groupArray[chain][0]][groupArray[chain][1]][groupArray[chain][2]] < 100: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.3f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][0]][groupArray[chain][1]][groupArray[chain][2]]) ) + elif self.__density[groupArray[chain][0]][groupArray[chain][1]][groupArray[chain][2]] < 1000: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.2f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][0]][groupArray[chain][1]][groupArray[chain][2]]) ) + elif self.__density[groupArray[chain][0]][groupArray[chain][1]][groupArray[chain][2]] < 10000: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.1f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][0]][groupArray[chain][1]][groupArray[chain][2]]) ) + else: + pdbFH.write ("ATOM %5i %4i DEN a 1 %8.3f%8.3f%8.3f 1.00%6.0f D\n" %( atomCount,chain,x, y, z, self.__density[groupArray[chain][0]][groupArray[chain][1]][groupArray[chain][2]]) ) + + atomCount += 1 +### chain should not be bigger than 99999 ... ehm shall be handeld in following versions .... *hust hust* + chain += 1 + pdbFH.write("TER\n") + pdbFH.write("END\n") + + pdbFH.close + print "saved "+filename diff --git a/lib/writepdb.pyc b/lib/writepdb.pyc new file mode 100644 index 0000000..61e4752 Binary files /dev/null and b/lib/writepdb.pyc differ diff --git a/pdb2dx.py b/pdb2dx.py new file mode 100755 index 0000000..1090d54 --- /dev/null +++ b/pdb2dx.py @@ -0,0 +1,82 @@ +#!/usr/bin/env python +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +### +# This script converts PDB files into OpenDX files, optionally a previously group can be set. Use this tool only for pdb files dxtuber previously created. +# +# v0.084 +# +settings are stored now in pdb files +# +filter by density + +import sys, os +from lib.tubedx import * +from lib.buildopendx import * +from lib.readpdb import * +from lib.writedx import * + + + +dxTuberConfFile = os.environ['HOME']+"/.dxTuber/dxTuber.conf" +conf = open (dxTuberConfFile,'r') +confLines = conf.readlines() +for line in confLines: + if re.match ('^#',line): + continue + if re.match('VERSION',line): +# print line + tmpArray = re.split('\s',line,1) + tmp = re.sub('\"','',tmpArray[1]) + tmp = re.sub('\n','',tmp) + version = re.sub('^\s','',tmp) +print "dxTuber v" + version + + +try: + sys.argv[1] #pdb file + sys.argv[2] #dx file +except IndexError: + print "\nThis script converts PDB files into OpenDX files, optionally a group ca be given as argument .\n\nusage:\npdb2dx (group)\n<> Mandatory parameters \n() Optional parameters\n\n" + sys.exit() +try: + sys.argv[3] + group = sys.argv[3] +except: + group = 0 + +print "Reading..." +readPDB = ReadPDB (sys.argv[1]) +if readPDB.readPDB() == 0: + dxObject = readPDB.getDxObject() +else: + print "Can not read pdb file" + sys.exit() +# only selected group gets densities > 0 +if group: + print "Extracting group id: "+ sys.argv[3] + voxelGroupGrid = readPDB.getTubeDxObject().getVoxelGroupGrid() #[x][y][z] = group + densityArray = dxObject.getDensity() #[x][y][z] = density + for xCount in range (len (densityArray)): + for yCount in range (len (densityArray[0])): + for zCount in range (len (densityArray[0][0])): + if sys.argv[3] != voxelGroupGrid[xCount][yCount][zCount]: + densityArray[xCount][yCount][zCount]=0 + dxObject.setDensity(densityArray) +saveToDx = WriteDx (dxObject) +saveToDx.write(sys.argv[2]) +print "Done" + + diff --git a/postgroupcavities.py b/postgroupcavities.py new file mode 100755 index 0000000..3ad563b --- /dev/null +++ b/postgroupcavities.py @@ -0,0 +1,164 @@ +#!/usr/bin/env python +# +# This file is part of the dxTuber package +# Copyright (C) 2010-2013 Martin Raunest +# +# 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 2 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, write to the Free Software Foundation, Inc., +# 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +# + + +# +# v0.23 +# +distcutoff parameter can now be set via cmd +# +# +# v0.2 +# + shell script to postgroup cavities into subcavities. +# +# +# + +import os +import subprocess +from lib.readpdb import * +from lib.writepdb import * + +from lib.buildopendx import * +from lib.tubedx import * + +from lib.postgroup import * + + + +######### +# Get dxTuber settings +######### + + + +dxTuberConfFile = os.environ['HOME']+"/.dxTuber/dxTuber.conf" +conf = open (dxTuberConfFile,'r') +confLines = conf.readlines() + +for line in confLines: + if re.match ('^#',line): + continue + if re.match('VERSION',line): + tmpArray = re.split('\s',line,1) + tmp = re.sub('\"','',tmpArray[1]) + tmp = re.sub('\n','',tmp) + version = re.sub('^\s','',tmp) + +#### +# Deal with arguments +#### + +###################################################################### +# If argparse is not installed global, get a static lib # +# and import it manually (static) # +# # +# http://packages.ubuntu.com/lucid/python-argparse # +# http://code.google.com/p/argparse/ # +# # +# maybe this helps for implementation ... # +# # +# sys.path.append('path_to_argparse_lib') +# from argparse import * +# # +###################################################################### + + +try: + import argparse + parser = argparse.ArgumentParser(description='''"--cubicsize" defines the number of voxels in each direction, e.g: 2 => 2 in +x and -x, +y and -y, +z and -z direction defining a 5x5x5 cube with 125 voxels + +"-t" sets the minimum amount of neighboars in the "--cubicsize" defined cube as a threshold for the core cavities''', + epilog='There is no spoon',formatter_class=argparse.RawTextHelpFormatter ) +except ImportError : + sys.path.append('/home/bit/raunest/programms/python/scripts/lib') ## if argparse is not installed global you can use it in a static way + from argparse import * + parser = ArgumentParser(description='"--cubicsize" defines the number of voxels in each direction, e.g: 2 => 2 in +x and -x, +y and -y, +z and -z direction defining a 5x5x5 cube with 125 voxels\n\n"-t" sets the minimum amount of neighboars in the "--cubicsize" defined cube as a threshold for the core cavities',epilog='There is no spoon',formatter_class=RawTextHelpFormatter ) + +parser.add_argument('--ipdb', default=False,help='PDB inputfile (only previously via dxTuber created pdb files are supported)') +parser.add_argument('--opdb', default=False, help='PDB outfile') +parser.add_argument('--cid', default='all', help='Cavity ID to postgroup, use "all" if all cavities should be processed (default="all")') +parser.add_argument('--cubicsize', default='2', help='edge length default 2 Voxels => Volume => 125 Voxels ... [INT] ') +parser.add_argument('-t','--threshold', default='100', help='Minimum amount of neighbors in --cubicsize to calc "core" cavities (default = 100)') +parser.add_argument('-d','--distcutoff', default='8', help='Cubic size for reprouping the non core voxels. This routine checks each core voxel per default in +8 and -8 in xyz for non core ISVs.') + +args = parser.parse_args() + +if not args.ipdb or not args.opdb: + sysRun = subprocess.Popen([sys.argv[0], '-h']) + sysRun.communicate() + print "Please specify input and output by '--ipdb' and '--opdb'\n" + sys.exit() + + + +readPDB = ReadPDB(args.ipdb) +readPDB.readPDB() +cavitiesTubeObject = readPDB.getTubeDxObject() +cavitiesDxObject = readPDB.getDxObject() + +### Densities contain now information about each voxels " neighborhood " +postGroup = PostGroup ( cavitiesTubeObject ) + +#oke testing now :) + + +# postGroupDensities = postGroup.calcNeighbors(cavityID = args.cid, edgeLength = int(args.cubicsize)) <-- works fine ! +# postGroupDxObject = BuildOpenDx (density = postGroupDensities, +# stepSize = cavitiesDxObject.getStepsize(), +# origin = cavitiesDxObject.getOrigin(), +# dimention = cavitiesDxObject.getDimention()) +# + + + + + + +''' +postGroupTubeDx = TubeDx ( + openDxObject = postGroupDxObject, + filename = cavitiesTubeObject.getFilename(), + proteinObject= cavitiesTubeObject.getProteinObject(), + waterObject= cavitiesTubeObject.getWaterObject(), + scanMethod = cavitiesTubeObject.getScanMethod(), + scanned=cavitiesTubeObject.getScanned(), + grouped='Yes', + groupes=None, + VoxelGroupGrid = None, + minDiameter = cavitiesTubeObject.getMinDiameter(), + protThreshold = cavitiesTubeObject.getProtThreshold(), + solventThreshold =cavitiesTubeObject.getSolventThreshold(), + version = version, + protFile =cavitiesTubeObject.getProtFile(), + solvFile = cavitiesTubeObject.getSolvFile(), + filterApplied = cavitiesTubeObject.getFilterApplied(), + ) +''' +## tested ! until here .... + +postGroupTubeDx = postGroup.postGroup_Neighbor(threshold =int(args.threshold), cavityID = args.cid, edgeLength = int(args.cubicsize) ,distCutoff =int(args.distcutoff)) + +writePDB = WritePDB (postGroupTubeDx) +#writePDB.write(filename= "2"+args.opdb) +writePDB.writeGroups(groupArray= postGroupTubeDx.getGroupes(), filename= args.opdb) + +print "Done" + +