| > |
restart: interface (warnlevel=0): with(plottools): with(StringTools): with(plots): with(Maplets[Elements]):
# +++++++++++++++++++
# MyMaplet creates a maplet window.
# Indentation is unnecessary, but hopefully helps to
# clarify syntax and the long list of parameters.
# +++++++++++++++++++
MyMaplet := Maplet
(
'onstartup'= RunWindow('W1'),
Window['W1']
(
'title'="Plot Special PDB Files",
[
[
"Special PDB File:", TextField['TF1'](),
Button
(
"File...",
Action(Evaluate('TF1' = "SelectPDB"))
)
],
[
"Initial Aminoacid:", TextField['TF2']('value'=1),
"Final Aminoacid:", TextField['TF3']('value'=1),
"Plot Tags:", CheckBox['TF4']('value'='true')
],
[TextBox['TB1']('height'=5) ],
Plotter['PL1'](),
[
Button
(
"Plot",
Action(Evaluate('PL1' = PlotPDB(TF1,TF2,TF3,TF4)))
),
Button
(
"List Aminoacids",
Action(Evaluate('TB1' = LoadAmino(TF1)))
),
Button
(
"End", Shutdown()
)
]
]
),
FileDialog['filPDBFile']
(
'onapprove' = Shutdown(['filPDBFile']),
'oncancel' = Shutdown()
)
):
# +++++++++++++++++++
# Procedure SelectPDB
# +++++++++++++++++++
SelectPDB := proc()
file := Maplets[Display]( Maplet( FileDialog['FD1'](
'onapprove' = Shutdown(['FD1']),
'oncancel' = Shutdown()
) ) ):
s := file[1]:
t := SubstituteAll(s, "\\", "//");
end proc:
# +++++++++++++++++++
# Procedure LoadAmino
# +++++++++++++++++++
LoadAmino := proc(sFile)
# ************************************
# Load the protein in a data structure
# ************************************
lProtein := readdata(sFile,[string, integer, string, string, string, string, float, float, float, float, float, string]):
# ******************************************************
# Store the atom range of every aminoacid in the protein
# ******************************************************
sCurAmn := lProtein[1][4]:
lCurAmn := lProtein[1][6]:
nCurAmn := 1:
lAminoacid := []:
n:=0:
for i from 1 to nops(lProtein) do
if i < nops(lProtein) then
sNxtAmn := lProtein[i][4]:
lNxtAmn := lProtein[i][6]:
if lCurAmn <> lNxtAmn then
n := n+1:
lAminoacid := [op(lAminoacid),[n, sCurAmn, nCurAmn, i-1]]:
lCurAmn := lNxtAmn:
sCurAmn := sNxtAmn:
nCurAmn := i:
end if:
else
lAminoacid := [op(lAminoacid),[n+1, sCurAmn, nCurAmn, i-1]]:
end if:
end do:
lAminoacid;
end proc:
# +++++++++++++++++
# Procedure PlotPDB
# +++++++++++++++++
PlotPDB := proc(sFile, nAmnIni, nAmnEnd, bShowTags)
# Avoid the most likely to happen error: nAmnIni > nAmnEnd
if nAmnIni > nAmnEnd then
newAmnEnd := nAmnIni:
else
newAmnEnd := nAmnEnd:
end if:
# ************************************
# Load the protein in a data structure
# ************************************
lProtein := readdata(sFile,[string, integer, string, string, string, string, float, float, float, float, float, string]):
# ******************************************************
# Store the atom range of every aminoacid in the protein
# ******************************************************
sCurAmn := lProtein[1][4]:
lCurAmn := lProtein[1][6]:
nCurAmn := 1:
lAminoacid := []:
for i from 1 to nops(lProtein) do
if i < nops(lProtein) then
sNxtAmn := lProtein[i][4]:
lNxtAmn := lProtein[i][6]:
if lCurAmn <> lNxtAmn then
lAminoacid := [op(lAminoacid),[sCurAmn, nCurAmn, i-1]]:
lCurAmn := lNxtAmn:
sCurAmn := sNxtAmn:
nCurAmn := i:
end if:
else
lAminoacid := [op(lAminoacid),[sCurAmn, nCurAmn, i-1]]:
end if:
end do:
# ********************************
# Initialize the lists for graphics
# ********************************
lCoordList := []:
lLabelList := []:
lCarbon := []:
lOxigen := []:
lNitrogen := []:
lSulfur := []:
lPolyAmn := []:
lAminoNameList := []:
# ******************************************************
# Build the polygons for the aminoacids
# ******************************************************
for j from nAmnIni to newAmnEnd do
i := lAminoacid[j][2]:
k := lAminoacid[j][3]:
PolySideChain := []:
# W is the first atom of the next aminoacid. Used for graphic continuity
W:=[lProtein[k+1][7],lProtein[k+1][8],lProtein[k+1][9]]:
# The first 4 atoms are always the same for every aminoacid in the sample PDB file
N := [lProtein[i][7],lProtein[i][8],lProtein[i][9]]:
CA:= [lProtein[i+1][7],lProtein[i+1][8],lProtein[i+1][9]]:
C := [lProtein[i+2][7],lProtein[i+2][8],lProtein[i+2][9]]:
OX:= [lProtein[i+3][7],lProtein[i+3][8],lProtein[i+3][9]]:
PolyBackbone := CURVES([N,CA,C,OX,C,W], THICKNESS(3));
# Find the geometric center of the aminoacid and put there the three-letter name
if k-i > 4 then
X := sum('lProtein[l][7]','l'=i+4..k)/(k-i-3):
Y := sum('lProtein[l][8]','l'=i+4..k)/(k-i-3):
# A little displacement in Z needed for Asparagine
Z := sum('lProtein[l][9]','l'=i+4..k)/(k-i-3) - 0.4:
else
# This is the case of the Glycine, which doesn't have a side chain
# and Alaninie, which just have one carbon as a side chain.
# For those two, we put the tag next to the alpha carbon.
X := lProtein[i+1][7] + 0.3:
Y := lProtein[i+1][8] + 0.3:
Z := lProtein[i+1][9] + 0.3:
end if:
AminoGeomCentr := [X, Y, Z, lAminoacid[j][1]]:
lAminoNameList:= [op(lAminoNameList),AminoGeomCentr]:
# ******************************************************
# Aliphatic aminoacids (1-4)
# ******************************************************
# 1) A = ALA = Alanine
if lProtein[i][4] = "ALA" then
CB:= [lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
PolySideChain := CURVES([CA,CB]);
end if:
# 2) V = VAL = Valine
if lProtein[i][4] = "VAL" then
CB:= [lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG1:=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
CG2:=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
PolySideChain := CURVES([CA,CB,CG2,CB,CG1]);
end if:
# 3) L = LEU = Leucine
if lProtein[i][4] = "LEU" then
CB:= [lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG:=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
CD1:=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
CD2:=[lProtein[i+7][7],lProtein[i+7][8],lProtein[i+7][9]]:
PolySideChain := CURVES([CA,CB,CG,CD1,CG,CD2]);
end if:
# 4) I = ILE = Isoleucine
if lProtein[i][4] = "ILE" then
CB:= [lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG1:=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
CG2:=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
CD1:=[lProtein[i+7][7],lProtein[i+7][8],lProtein[i+7][9]]:
PolySideChain := CURVES([CA,CB,CG2,CB,CG1,CD1]);
end if:
# ******************************************************
# Hydrophobic-aromatic aminoacids (5-7)
# ******************************************************
# 5) F = PHE = Phenylanaline
if lProtein[i][4] = "PHE" then
CB:= [lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG:=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
CD1:=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
CD2:=[lProtein[i+7][7],lProtein[i+7][8],lProtein[i+7][9]]:
CE1:=[lProtein[i+8][7],lProtein[i+8][8],lProtein[i+8][9]]:
CE2:=[lProtein[i+9][7],lProtein[i+9][8],lProtein[i+9][9]]:
CZ:=[lProtein[i+10][7],lProtein[i+10][8],lProtein[i+10][9]]:
PolySideChain := CURVES([CA,CB,CG,CD1,CE1,CZ,CE2,CD2,CG]);
end if:
# 6) Y = TYR = Tyrosine
if lProtein[i][4] = "TYR" then
CB:= [lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG:=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
CD1:=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
CD2:=[lProtein[i+7][7],lProtein[i+7][8],lProtein[i+7][9]]:
CE1:=[lProtein[i+8][7],lProtein[i+8][8],lProtein[i+8][9]]:
CE2:=[lProtein[i+9][7],lProtein[i+9][8],lProtein[i+9][9]]:
CZ:=[lProtein[i+10][7],lProtein[i+10][8],lProtein[i+10][9]]:
OH:=[lProtein[i+11][7],lProtein[i+11][8],lProtein[i+11][9]]:
PolySideChain := CURVES([CA,CB,CG,CD1,CE1,CZ,OH,CZ,CE2,CD2,CG]);
end if:
# 7) W = TRP = Tryptophan
if lProtein[i][4] = "TRP" then
CB:= [lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG:=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
CD1:=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
CD2:=[lProtein[i+7][7],lProtein[i+7][8],lProtein[i+7][9]]:
NE1:=[lProtein[i+8][7],lProtein[i+8][8],lProtein[i+8][9]]:
CE2:=[lProtein[i+9][7],lProtein[i+9][8],lProtein[i+9][9]]:
CE3:=[lProtein[i+10][7],lProtein[i+10][8],lProtein[i+10][9]]:
CZ2:=[lProtein[i+11][7],lProtein[i+11][8],lProtein[i+11][9]]:
CZ3:=[lProtein[i+12][7],lProtein[i+12][8],lProtein[i+12][9]]:
CH2:=[lProtein[i+13][7],lProtein[i+13][8],lProtein[i+13][9]]:
PolySideChain := CURVES([CA,CB,CG,CD1,NE1,CE2,CZ2,CH2,CZ3,CE3,CD2,CE2,CD2,CG]);
end if:
# ******************************************************
# Neutral-polar aminoacids (8-13)
# ******************************************************
# 8) S = SER = Serine
if lProtein[i][4] = "SER" then
CB:= [lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
OG:=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
PolySideChain := CURVES([CA,CB,OG]);
end if:
# 9) T = THR = Threionine
if lProtein[i][4] = "THR" then
CB:= [lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
OG1:=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
CG2:=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
PolySideChain := CURVES([CA,CB,CG2,CB,OG1]);
end if:
# 10) C = CYS = Cysteine
if lProtein[i][4] = "CYS" then
CB:= [lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
SG:=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
PolySideChain := CURVES([CA,CB,SG]);
end if:
# 11) M = MET = Methionine
if lProtein[i][4] = "MET" then
CB:= [lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG:=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
SD:=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
CE:=[lProtein[i+7][7],lProtein[i+7][8],lProtein[i+7][9]]:
PolySideChain := CURVES([CA,CB,CG,SD,CE]);
end if:
# 12) N = ASN = Asparagine
if lProtein[i][4] = "ASN" then
CB :=[lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG :=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
OD1:=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
ND2:=[lProtein[i+7][7],lProtein[i+7][8],lProtein[i+7][9]]:
PolySideChain := CURVES([CA,CB,CG,OD1,CG,ND2]);
end if:
# 13) Q = GLN = Glutamine
if lProtein[i][4] = "GLN" then
CB := [lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG :=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
CD :=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
OE1:=[lProtein[i+7][7],lProtein[i+7][8],lProtein[i+7][9]]:
NE2:=[lProtein[i+8][7],lProtein[i+8][8],lProtein[i+8][9]]:
PolySideChain := CURVES([CA,CB,CG,CD,OE1,CD,NE2]);
end if:
# ******************************************************
# Aliphatic aminoacids (14-15)
# ******************************************************
# 14) D = ASP = Aspartate
if lProtein[i][4] = "ASP" then
CB :=[lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG :=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
OD1:=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
OD2:=[lProtein[i+7][7],lProtein[i+7][8],lProtein[i+7][9]]:
PolySideChain := CURVES([CA,CB,CG,OD1,CG,OD2]);
end if:
# 15) E = GLU = Glutamate
if lProtein[i][4] = "GLU" then
CB :=[lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG :=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
CD :=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
OE1:=[lProtein[i+7][7],lProtein[i+7][8],lProtein[i+7][9]]:
OE2:=[lProtein[i+8][7],lProtein[i+8][8],lProtein[i+8][9]]:
PolySideChain := CURVES([CA,CB,CG,CD,OE1,CD,OE2]);
end if:
# ******************************************************
# Basic aminoacids (16-18)
# ******************************************************
# 16) H = HIS = Histidine
if lProtein[i][4] = "HIS" then
CB :=[lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG :=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
ND1:=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
CD2:=[lProtein[i+7][7],lProtein[i+7][8],lProtein[i+7][9]]:
CE1:=[lProtein[i+8][7],lProtein[i+8][8],lProtein[i+8][9]]:
NE2:=[lProtein[i+9][7],lProtein[i+9][8],lProtein[i+9][9]]:
PolySideChain := CURVES([CA,CB,CG,ND1,CE1,NE2,CD2,CG]);
end if:
# 17) K = LYS = Lysine
if lProtein[i][4] = "LYS" then
CB:=[lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG:=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
CD:=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
CE:=[lProtein[i+7][7],lProtein[i+7][8],lProtein[i+7][9]]:
NZ:=[lProtein[i+8][7],lProtein[i+8][8],lProtein[i+8][9]]:
PolySideChain := CURVES([CA,CB,CG,CD,CE,NZ]);
end if:
# 18) R = ARG = Arginine
if lProtein[i][4] = "ARG" then
CB :=[lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG :=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
CD :=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
NE :=[lProtein[i+7][7],lProtein[i+7][8],lProtein[i+7][9]]:
CZ :=[lProtein[i+8][7],lProtein[i+8][8],lProtein[i+8][9]]:
NH1:=[lProtein[i+9][7],lProtein[i+9][8],lProtein[i+9][9]]:
NH2:=[lProtein[i+10][7],lProtein[i+10][8],lProtein[i+10][9]]:
PolySideChain := CURVES([CA,CB,CG,CD,NE,CZ,NH1,CZ,NH2]);
end if:
# ******************************************************
# Conformationally important aminoacids (19-20)
# ******************************************************
# 19) G = GLY = Glycine
if lProtein[i][4] = "GLY" then
# Glycine doesn't have a side chain
end if:
# 20) P = PRO = Proline
if lProtein[i][4] = "PRO" then
CB:=[lProtein[i+4][7],lProtein[i+4][8],lProtein[i+4][9]]:
CG:=[lProtein[i+5][7],lProtein[i+5][8],lProtein[i+5][9]]:
CD:=[lProtein[i+6][7],lProtein[i+6][8],lProtein[i+6][9]]:
PolySideChain := CURVES([CA,CB,CG,CD,N]);
end if:
# Add sidechain, backbone and aminoacid name to the list of itme sto be plotted
if nops(PolySideChain) > 0 then
lPolyAmn := [op(lPolyAmn),PolySideChain]:
end if:
lPolyAmn := [op(lPolyAmn),PolyBackbone]:
end do:
# ******************************************************************************
# Draw the atom's number and solids of different colors depending upon the atoms
# ******************************************************************************
for i from lAminoacid[nAmnIni][2] to lAminoacid[newAmnEnd][3] do
lCoord := [lProtein[i][7], lProtein[i][8], lProtein[i][9]]:
lLabel := [lProtein[i][7], lProtein[i][8], lProtein[i][9]+0.4, lProtein[i][3]]:
if lProtein[i][12] = "C" then
lCarbon := [op(lCarbon),lCoord]:
end if:
if lProtein[i][12] = "O" then
lOxigen := [op(lOxigen),lCoord]:
end if:
if lProtein[i][12] = "N" then
lNitrogen := [op(lNitrogen),lCoord]:
end if:
if lProtein[i][12] = "S" then
lSulfur := [op(lSulfur),lCoord]:
end if:
lCoordList := [op(lCoordList),lCoord]:
lLabelList := [op(lLabelList),lLabel]:
end do:
# ******************************************************
# Load all the objects into the PLO3D structure
# ******************************************************
polCarbon:=polyhedraplot(lCarbon,polytype=hexahedron,polyscale =0.2,lightmodel=light4,shading=XYZ,style=PATCHNOGRID,color=gray):
polOxygen:=polyhedraplot(lOxigen,polytype=hexahedron,polyscale =0.2,lightmodel=light4,shading=XYZ,style=PATCHNOGRID,color=red):
polNitrogen:=polyhedraplot(lNitrogen,polytype=hexahedron,polyscale =0.2,lightmodel=light4,shading=XYZ,style=PATCHNOGRID,color=blue):
AtomLabel := textplot3d(lLabelList,color=black):
AminoLabel := textplot3d(lAminoNameList,color=red):
lSideChains := PolySideChains(lAminoacid):
if nops(lSulfur) > 0 then
polSulfur:=polyhedraplot(lSulfur, polytype=hexahedron,
polyscale=0.3, lightmodel=light4, shading=XYZ,
style=PATCHNOGRID, color=yellow):
if bShowTags = 'true' then
display3d(lPolyAmn, AtomLabel, AminoLabel, polCarbon,
polOxygen, polNitrogen, polSulfur, color=black);
else
display3d(lPolyAmn, polCarbon, polOxygen, polNitrogen, polSulfur, color=black);
end if:
else
if bShowTags = 'true' then
display3d(lPolyAmn, AtomLabel, AminoLabel, polCarbon,
polOxygen, polNitrogen, color=black);
else
display3d(lPolyAmn, polCarbon, polOxygen, polNitrogen, color=black);
end if:
end if:
end proc:
Maplets[Display](MyMaplet); |
