Hi Thomas,
You're right, and using an svd on the whole ring is also more correct.
But the vector stuff is simple and handy to have lying around, and it
serves also to give some insight in what's going on :)
numpy.linalg.svd is magic if you don't know how a cross product or
normalization works... So, now to make it a function taking two ring
bearing residues as arguments:
def ring_angle(selection1, selection2):
a = plane_normal(selection1 + " and not (e. h or n. c,n,ca,cb,o)")
b = plane_normal(selection2 + " and not (e. h or n. c,n,ca,cb,o)")
return math.degrees(cpv.get_angle(a, b))
cmd.extend("ring_angle",ring_angle)
Ramiro, do we need to write the methods section for the manuscript too? :p
Cheerio :)
Tsjerk
On Thu, Sep 23, 2010 at 1:33 PM, Thomas Holder
<spel...@users.sourceforge.net> wrote:
> no need to implement common linear algebra functions, there is the
> chempy.cpv module shipped with pymol (and there is numpy as well).
>
> Ramiro, I recently did something similar, just adjust the residue
> selection in the code below (requires numpy):
>
> python
> from chempy import cpv
> import numpy, math
> def plane_normal(selection):
> stored.x = list()
> cmd.iterate_state(-1, selection, 'stored.x.append([x,y,z])')
> x = numpy.array(stored.x)
> U,s,Vh = numpy.linalg.svd(x - x.mean(0))
> return cpv.normalize(Vh[2])
> dir1 = plane_normal('A/37/CG+CD1+CE1+CZ+CE2+CD2')
> dir2 = plane_normal('A/41/CG+CD1+CE1+CZ+CE2+CD2')
> print 'Angle in degrees:', math.degrees(cpv.get_angle(dir1, dir2))
> python end
>
> Cheers,
> Thomas
>
> On Thu, 2010-09-23 at 13:02 +0200, Tsjerk Wassenaar wrote:
>> Hi Ramiro,
>>
>> A bit of linear algebra wouldn't hurt... :p
>> In python:
>>
>> def vsub(a,b): return a[0]-b[0], a[1]-b[1], a[2]-b[2]
>>
>> def dot(a,b): return a[0]*b[0]+a[1]*b[1]+a[2]*b[2]
>>
>> def svmul(s,a): return s*a[0], s*a[1], s*a[2]
>>
>> def normalize(a): return svmul(1/math.sqrt(dot(a,a)),a)
>>
>> def cross(a,b): return a[1]*b[2]-a[2]*b[1], a[2]*b[0]-a[0]*b[2],
>> a[0]*b[1]-a[1]*b[0]
>>
>> a = cmd.get_model('r. phe and i. ##RESIDUE1## and n. cg,ce1,ce2').atom
>> a = [ i.coord for i in a ]
>> b = cmd.get_model('r. phe and i. ##RESIDUE2## and n. cg,ce1,ce2').atom
>> b = [ i.coord for i in b ]
>>
>> na = normalize(cross(vsub(a[1],a[0]),vsub(a[2],a[0])))
>> nb = normalize(cross(vsub(b[1],b[0]),vsub(b[2],b[0])))
>> angle = math.acos(dot(na,nb))
>>
>> print angle
>>
>> ###
>>
>> Haven't tested it, and there may be more efficient ways of getting the
>> coordinates. If you run into problems like this more often, it's
>> likely that you should pick up on algebra and programming... :)
>>
>> Have fun,
>>
>> Tsjerk
>>
>>
>> On Thu, Sep 23, 2010 at 12:19 PM, Ramiro Téllez Sanz
>> <urcind...@gmail.com> wrote:
>> > Thanks for your kind help, Tsjerk.
>> >
>> >> Hi Ramiro,
>> >>
>> >> Assuming your rings are nicely planar, and representing the ring as:
>> >>
>> >> 1-2-3
>> >> | |
>> >> 6-5-4
>> >>
>> >> you can get the plane normal vector as the vector cross product from
>> >> (3)-(1) and (5)-(1).
>> >
>> > OK. But I just started to use pymol. Which are the commands to do so?
>> > I know how to get the coordinates of a selected atom, but need the pymol
>> > commands to treat the data:
>> > a) How to create the vectors from 1->3 and from 1->5
>> > b) How to treat the vectors to perform the vector cross product
>> >
>> >> Doing so for both rings gives you the two normal vectors. The angle
>> >> then follows from the dot product of the (normalized) normal vectors:
>> >>
>> >> angle = acos(n1 . n2)
>> >
>> > Again, I would need the commands to:
>> > c) Normalize the vectors (how to set their modules = 1)
>> >
>> > I also guess n1 and n2 represent the normalized vectors, don't they? So
>> > this command is very clear :)
>> >
>> >> It becomes a bit more elaborate if the planes are not planar :)
>> >>
>> >> Hope it helps,
>> >>
>> >> Tsjerk
>> >
>> > Again, thanks very much in advance for your kind help.
>> >
>> >> On Thu, Sep 23, 2010 at 10:53 AM, Ramiro Téllez Sanz
>> >> <urcind...@gmail.com> wrote:
>> >>>
>> >>> Hi everyone and thanks for reading this!
>> >>>
>> >>> I am interested in measuring the angle between aromatic ring planes.
>> >>> Is there any easy way/script to do it?
>> >>>
>> >>> One way that came to my mind is creating a pseudoatom representing the
>> >>> centroid for each ring (I already know how to do that), then drawing two
>> >>> lines perpendicularly to the planes from both centroids, and finally
>> >>> measuring the angle between the lines. Will that be possible? How could
>> >>> this be done?
>> >>>
>> >>> Is there any other way? I'm completely clueless. Any help will be
>> >>> greatly appreciated.
>> >>>
>> >>> Thanks in advance,
>> >>>
>> >>> Ramiro Tellez Sanz
>> >>> Dept. Physical Chemistry
>> >>> University of Almeria
>> >>> Spain
>>
>>
>> --
>> Tsjerk A. Wassenaar, Ph.D.
>>
>> post-doctoral researcher
>> Molecular Dynamics Group
>> Groningen Institute for Biomolecular Research and Biotechnology /
>> University of Groningen
>> The Netherlands
> --
> Thomas Holder
> Group of Steffen Schmidt
> Department of Biochemistry
> MPI for Developmental Biology
> Spemannstr. 35
> D-72076 Tübingen
>
>
>
>
> ------------------------------------------------------------------------------
> Start uncovering the many advantages of virtual appliances
> and start using them to simplify application deployment and
> accelerate your shift to cloud computing.
> http://p.sf.net/sfu/novell-sfdev2dev
> _______________________________________________
> PyMOL-users mailing list (PyMOL-users@lists.sourceforge.net)
> Info Page: https://lists.sourceforge.net/lists/listinfo/pymol-users
> Archives: http://www.mail-archive.com/pymol-users@lists.sourceforge.net
--
Tsjerk A. Wassenaar, Ph.D.
post-doctoral researcher
Molecular Dynamics Group
* Groningen Institute for Biomolecular Research and Biotechnology
* Zernike Institute for Advanced Materials
University of Groningen
The Netherlands
------------------------------------------------------------------------------
Start uncovering the many advantages of virtual appliances
and start using them to simplify application deployment and
accelerate your shift to cloud computing.
http://p.sf.net/sfu/novell-sfdev2dev
_______________________________________________
PyMOL-users mailing list (PyMOL-users@lists.sourceforge.net)
Info Page: https://lists.sourceforge.net/lists/listinfo/pymol-users
Archives: http://www.mail-archive.com/pymol-users@lists.sourceforge.net
