Avogadro has a plugin architecture that allows you to drop a Python script into a particular folder to add additional functionality, whether it be a new menu item, a new display method or whatever. I've previously used this to write
an extension to read comp chem log files using cclib. This time I need to add a new display method (or Engine in the lingo of Avogadro) that will depict a pharmacophore.
So what's
Pharmer? Well, until the release of Pharao by
Silicos last year there were no open source pharmacophore programs. And the recent release of
Pharmer by David Koes bumps the number up to two. Pharmer focuses on pharmacophore search - ultra-fast pharmacophore search, that is. It does this by precomputing an index that allows very efficient searching. The common denominator of Pharmer and Pharao is not only "Phar", but also that both use the
Open Babel library.
The file format used by Pharmer is simply JSON, which is easily edited by a text editor. Here's an example that was created using the
ZINC Pharmer interface:
HSP90.txt. If you install Avogadro 1.0.1 with Python support on Windows, also have Python 2.6 installed in C:\Python26, and copy the Python script shown below to C:\Program Files (x86)\Avogadro\bin\engineScripts, then you can add the Pharmacophore engine (click "Add" in the "Display Types" panel; note: you can add more than one Pharmacophore engine if you want to display more than one pharmacophore at the same time). Once added, click on Settings and Browse to open the Pharmer JSON file. If nothing appears on the screen, try rotating the display slightly.
Note: The image above was created by joining the output of the following Pharmer command using "obabel --join" and depicting with Avogadro:
pharmer dbsearch -dbdir db_712_actives -in HSP90.txt -out results.sdf -reduceconfs 1 -sort-rmsd
pharmacophore.py:
import os
from PyQt4.Qt import *
import Avogadro as avo
from numpy import array
# Get JSON from my own Python installation
import sys
sys.path.append("C:\\Python26\\Lib\\site-packages")
sys.path.append("C:\\Python26\\lib")
import json
colorlookup = {
'Hydrophobic': [0, 139, 69],
'HydrogenAcceptor': [255, 165, 0],
'HydrogenDonor': [211, 211, 211],
'Aromatic': [205, 0, 205],
'PositiveIon': [0, 0, 238],
'NegativeIon': [255, 0, 0],
}
class Engine(QObject):
# declare the changed() signal
__pyqtSignals__ = ("changed()",)
# constructor
def __init__(self):
QObject.__init__(self)
self.widget = None
self.textbox = None
self.filename = ""
self.pharma = []
def name(self):
return "Pharmacophore"
def flags(self):
return avo.EngineFlags.NoFlags
def debug(self, text):
return
filename = os.path.join("C:\\", "Users", "Noel", "Desktop", "tmp.txt")
output = open(filename, "a")
print >> output, text
output.close()
@pyqtSignature("")
def browse(self):
self.filename = str(QFileDialog.getOpenFileName())
if not self.filename:
return
self.textbox.setText(self.filename)
with open(self.filename) as f:
query = json.load(f)
self.pharma = [x for x in query['points'] if x['enabled']]
if len(self.pharma) > 0:
for p in self.pharma:
x += p['x']
y += p['y']
z += p['z']
x /= len(self.pharma)
y /= len(self.pharma)
z /= len(self.pharma)
for p in self.pharma:
p['x'] -= x
p['y'] -= y
p['z'] -= z
self.debug("Pharma: %s" % self.pharma)
self.emit(SIGNAL("changed()"))
@pyqtSignature("")
def clear(self):
self.filename = ""
self.textbox.setText(self.filename)
self.pharma = []
self.emit(SIGNAL("changed()"))
def settingsWidget(self):
self.widget = QWidget()
layout = QVBoxLayout(self.widget)
self.textbox = QLineEdit(self.filename)
self.textbox.setReadOnly(True)
layout.addWidget(self.textbox)
buttons = QHBoxLayout(self.widget)
layout.addLayout(buttons)
self.browsebtn = QPushButton("Browse")
self.clearbtn = QPushButton("Clear")
buttons.addWidget(self.browsebtn)
buttons.addWidget(self.clearbtn)
QObject.connect(self.browsebtn, SIGNAL("clicked()"), self, SLOT("browse()"))
QObject.connect(self.clearbtn, SIGNAL("clicked()"), self, SLOT("clear()"))
return self.widget
def renderOpaque(self, pd):
self.debug("Len of pharma: %d" % len(self.pharma))
if len(self.pharma) == 0:
return
# Painter
painter = pd.painter
# Molecule
molecule = pd.molecule
# Color
color = pd.colorMap
for point in self.pharma:
if point['name'] in colorlookup:
colvals = [float(x) / 255. for x in colorlookup[point['name']]] + [1.0]
color.setFromRgba(*colvals)
else:
color.setFromRgba(0.3, 0.6, 1.0, 1.0)
painter.setColor(color)
r = point['radius']
self.debug("Radius: %.1f" % r)
x, y, z = point['x'], point['y'], point['z']
painter.drawSphere(array([x, y, z]), r)