path: root/cad/src/protein/commands/InsertPeptide/PeptideGenerator.py

# Copyright 2008 Nanorex, Inc.  See LICENSE file for details.
"""
PeptideGenerator.py

Peptide Generator can generate short polypeptide chains of arbitrarily choosen
sequence and secondary structure.

@author: Piotr
@version: $Id$
@copyright: 2008 Nanorex, Inc.  See LICENSE file for details.

@see http://www.nanoengineer-1.net/mediawiki/index.php?title=Peptide_generator_dialog
for notes about what's going on here.

History:

Ninad 2008-07-24: Refactoring / cleanup to port PeptideGenerator to the 
                 EditCommand API. (see InsertPeptide_EditCommand)
"""

import foundation.env as env

from geometry.InternalCoordinatesToCartesian import InternalCoordinatesToCartesian

from model.chem import Atom
from model.chunk import Chunk
from model.bond_constants import V_DOUBLE, V_AROMATIC
from operations.bonds_from_atoms import inferBonds

from protein.model.Protein import Protein 
from protein.model.Residue import Residue
from protein.model.Residue import SS_HELIX, SS_STRAND, SS_COIL, AA_3_TO_1

from Numeric import zeros, sqrt, pi, sin, cos, Float
from geometry.VQT import Q, V, norm, vlen, cross, angleBetween

from utilities.debug import print_compact_stack

# Internal coordinate sets for amino acids
# Converted from AMBER all_amino94.in file
# Fixed a few issues: conformation of phenylalanine, connectivity of serine.
# Beware of proline - it needs to be handled in a different way.

ALA_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp3",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB1", "H", "",       7,    5,    3,    1.090,  109.500,   60.000 ),
    (   9, "HB2", "H", "",       7,    5,    3,    1.090,  109.500,  180.000 ),
    (  10, "HB3", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (  11, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  12, "O  ", "O", "sp2",   11,    5,    3,    1.229,  120.500,    0.000 ),
]

GLY_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA2", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "HA3", "H", "",       5,    3,    2,    1.090,  109.500,   60.000 ),
    (   8, "C  ", "C", "sp2",    5,    3,    2,    1.522,  110.400,  180.000 ),
    (   9, "O  ", "O", "sp2",    8,    5,    3,    1.229,  120.500,    0.000 )
]

SER_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB2", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "HB3", "H", "",       7,    5,    3,    1.090,  109.500,   60.000 ),
    (  10, "OG ", "O", "sp3",    7,    5,    3,    1.430,  109.470,  180.000 ),
    (  11, "HG ", "H", "",      10,    7,    5,    0.960,  109.470,  180.000 ),
    (  12, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  13, "O  ", "O", "sp2",   12,    5,    3,    1.229,  120.500,    0.000 ),
]

PHE_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB2", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "HB3", "H", "",       7,    5,    3,    1.090,  109.500,   60.000 ),
    (  10, "CG ", "C", "sp2a",    7,    5,    3,    1.510,  115.000,  180.000 ),
    (  11, "CD1", "C", "sp2a",   10,    7,    5,    1.400,  120.000,  180.000 ),
    (  12, "HD1", "H", "",      11,   10,    7,    1.090,  120.000,    0.000 ),
    (  13, "CE1", "C", "sp2a",   11,   10,    7,    1.400,  120.000,  180.000 ),
    (  14, "HE1", "H", "",      13,   11,   10,    1.090,  120.000,  180.000 ),
    (  15, "CZ ", "C", "sp2a",   13,   11,   10,    1.400,  120.000,    0.000 ),
    (  16, "HZ ", "H", "",      15,   13,   11,    1.090,  120.000,  180.000 ),
    (  17, "CE2", "C", "sp2a",   15,   13,   11,    1.400,  120.000,    0.000 ),
    (  18, "HE2", "H", "",      17,   15,   13,    1.090,  120.000,  180.000 ),
    (  19, "CD2", "C", "sp2a",   17,   15,   13,    1.400,  120.000,    0.000 ),
    (  20, "HD2", "H", "",      19,   17,   15,    1.090,  120.000,  180.000 ),
    (  21, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  22, "O  ", "O", "sp2",   21,    5,    3,    1.229,  120.500,    0.000 ),
]

GLU_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB2", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "HB3", "H", "",       7,    5,    3,    1.090,  109.500,   60.000 ),
    (  10, "CG ", "C", "sp3",    7,    5,    3,    1.510,  109.470,  180.000 ),
    (  11, "HG2", "H", "",      10,    7,    5,    1.090,  109.500,  300.000 ),
    (  12, "HG3", "H", "",      10,    7,    5,    1.090,  109.500,   60.000 ),
    (  13, "CD ", "C", "sp2",   10,    7,    5,    1.527,  109.470,  180.000 ),
    (  14, "OE1", "O", "sp2",   13,   10,    7,    1.260,  117.200,   90.000 ),
    (  15, "OE2", "O", "sp3",   13,   10,    7,    1.260,  117.200,  270.000 ),
    (  16, "HE2", "H", "",      15,   13,   10,    0.960,  109.500,  180.000 ),
    (  17, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  18, "O  ", "O", "sp2",   17,    5,    3,    1.229,  120.500,    0.000 ),
]

PRO_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.337,  117.000,  180.000 ),
    (   4, "CD ", "C", "sp3",    3,    2,    1,    1.458,  126.100,  356.100 ),
    (   5, "HD2", "H", "",       4,    3,    2,    1.090,  109.500,   80.000 ),
    (   6, "HD3", "H", "",       4,    3,    2,    1.090,  109.500,  320.000 ),
    (   7, "CG ", "C", "sp3",    4,    3,    2,    1.500,  103.200,  200.100 ),
    (   8, "HG2", "H", "",       7,    4,    3,    1.090,  109.500,  218.000 ),
    (   9, "HG3", "H", "",       7,    4,    3,    1.090,  109.500,   98.000 ),
    (  10, "CB ", "C", "sp3",    7,    4,    3,    1.510,  106.000,  338.300 ),
    (  11, "HB2", "H", "",      10,    7,    4,    1.090,  109.500,  256.300 ),
    (  12, "HB3", "H", "",      10,    7,    4,    1.090,  109.500,  136.300 ),
    (  13, "CA ", "C", "sp3",    3,    2,    1,    1.451,  120.600,  175.200 ),
    (  14, "HA ", "H", "",      13,    3,    2,    1.090,  109.500,   60.000 ),
    (  15, "C  ", "C", "sp2",   13,    3,    2,    1.522,  109.500,  300.000 ),
    (  16, "O  ", "O", "sp2",   15,   13,    3,    1.229,  120.500,    0.000 ),
]

CYS_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB2", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "HB3", "H", "",       7,    5,    3,    1.090,  109.500,   60.000 ),
    (  10, "SG ", "S", "sp3",    7,    5,    3,    1.810,  116.000,  180.000 ),
    (  11, "HG ", "H", "",      10,    7,    5,    1.330,   96.000,  180.000 ),
    (  12, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  13, "O  ", "O", "sp2",   12,    5,    3,    1.229,  120.500,    0.000 ),
]

MET_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB2", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "HB3", "H", "",       7,    5,    3,    1.090,  109.500,   60.000 ),
    (  10, "CG ", "C", "sp3",    7,    5,    3,    1.525,  109.470,  180.000 ),
    (  11, "HG2", "H", "",      10,    7,    5,    1.090,  109.500,  300.000 ),
    (  12, "HG3", "H", "",      10,    7,    5,    1.090,  109.500,   60.000 ),
    (  13, "SD ", "S", "sp3",   10,    7,    5,    1.810,  110.000,  180.000 ),
    (  14, "CE ", "C", "sp3",   13,   10,    7,    1.780,  100.000,  180.000 ),
    (  15, "HE1", "H",  "",     14,   13,   10,    1.090,  109.500,   60.000 ),
    (  16, "HE2", "H",  "",     14,   13,   10,    1.090,  109.500,  180.000 ),
    (  17, "HE3", "H", "",      14,   13,   10,    1.090,  109.500,  300.000 ),
    (  18, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  19, "O  ", "O", "sp2",   18,    5,    3,    1.229,  120.500,    0.000 ),
]

THR_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB ", "H", "",       7,    5,    3,    1.090,  109.500,  180.000 ),
    (   9, "CG2", "C", "sp3",    7,    5,    3,    1.525,  109.470,  300.000 ),
    (  10, "HG2", "H", "",       9,    7,    5,    1.090,  109.500,   60.000 ),
    (  11, "HG2", "H", "",       9,    7,    5,    1.090,  109.500,  180.000 ),
    (  12, "HG2", "H", "",       9,    7,    5,    1.090,  109.500,  300.000 ),
    (  13, "OG1", "O", "sp3",    7,    5,    3,    1.430,  109.470,   60.000 ),
    (  14, "HG1", "H", "",      13,    7,    5,    0.960,  109.470,  180.000 ),
    (  15, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  16, "O  ", "O", "sp2",   15,    5,    3,    1.229,  120.500,    0.000 ),
]

LEU_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB2", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "HB3", "H", "",       7,    5,    3,    1.090,  109.500,   60.000 ),
    (  10, "CG ", "C", "sp3",    7,    5,    3,    1.525,  109.470,  180.000 ),
    (  11, "HG ", "H", "",      10,    7,    5,    1.090,  109.500,  300.000 ),
    (  12, "CD1", "C", "sp3",   10,    7,    5,    1.525,  109.470,   60.000 ),
    (  13, "HD1", "H", "",      12,   10,    7,    1.090,  109.500,   60.000 ),
    (  14, "HD1", "H", "",      12,   10,    7,    1.090,  109.500,  180.000 ),
    (  15, "HD1", "H", "",      12,   10,    7,    1.090,  109.500,  300.000 ),
    (  16, "CD2", "C", "sp3",   10,    7,    5,    1.525,  109.470,  180.000 ),
    (  17, "HD2", "H",  "",     16,   10,    7,    1.090,  109.500,   60.000 ),
    (  18, "HD2", "H",  "",     16,   10,    7,    1.090,  109.500,  180.000 ),
    (  19, "HD2", "H",  "",     16,   10,    7,    1.090,  109.500,  300.000 ),
    (  20, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  21, "O  ", "O", "sp2",   20,    5,    3,    1.229,  120.500,    0.000 ),
]


ILE_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  109.470,   60.000 ),
    (   8, "HB ", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "CG2", "C", "sp3",    7,    5,    3,    1.525,  109.470,   60.000 ),
    (  10, "HG2", "H", "",       9,    7,    5,    1.090,  109.500,   60.000 ),
    (  11, "HG2", "H", "",       9,    7,    5,    1.090,  109.500,  180.000 ),
    (  12, "HG2", "H", "",       9,    7,    5,    1.090,  109.500,  300.000 ),
    (  13, "CG1", "C", "sp3",    7,    5,    3,    1.525,  109.470,  180.000 ),
    (  14, "HG1", "H", "",      13,    7,    5,    1.090,  109.500,  300.000 ),
    (  15, "HG1", "H", "",      13,    7,    5,    1.090,  109.500,   60.000 ),
    (  16, "CD1", "C", "sp3",   13,    7,    5,    1.525,  109.470,  180.000 ),
    (  17, "HD1", "H",  "",     16,   13,    7,    1.090,  109.500,   60.000 ),
    (  18, "HD1", "H",  "",     16,   13,    7,    1.090,  109.500,  180.000 ),
    (  19, "HD1", "H",  "",     16,   13,    7,    1.090,  109.500,  300.000 ),
    (  20, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  21, "O  ", "O", "sp2",   20,    5,    3,    1.229,  120.500,    0.000 ),
]

VAL_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB ", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "CG1", "C", "sp3",    7,    5,    3,    1.525,  109.470,   60.000 ),
    (  10, "HG1", "H", "",       9,    7,    5,    1.090,  109.500,   60.000 ),
    (  11, "HG1", "H", "",       9,    7,    5,    1.090,  109.500,  180.000 ),
    (  12, "HG1", "H", "",       9,    7,    5,    1.090,  109.500,  300.000 ),
    (  13, "CG2", "C", "sp3",    7,    5,    3,    1.525,  109.470,  180.000 ),
    (  14, "HG2", "H", "",      13,    7,    5,    1.090,  109.500,   60.000 ),
    (  15, "HG2", "H", "",      13,    7,    5,    1.090,  109.500,  180.000 ),
    (  16, "HG2", "H", "",      13,    7,    5,    1.090,  109.500,  300.000 ),
    (  17, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  18, "O  ", "O", "sp2",   17,    5,    3,    1.229,  120.500,    0.000 ),
]

TRP_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB2", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "HB3", "H", "",       7,    5,    3,    1.090,  109.500,   60.000 ),
    (  10, "CG ", "C", "sp2",    7,    5,    3,    1.510,  115.000,  180.000 ),
    (  11, "CD1", "C", "sp2",   10,    7,    5,    1.340,  127.000,  180.000 ),
    (  12, "HD1", "H", "",      11,   10,    7,    1.090,  120.000,    0.000 ),
    (  13, "NE1", "N", "sp3",   11,   10,    7,    1.430,  107.000,  180.000 ),
    (  14, "HE1", "H", "",      13,   11,   10,    1.010,  125.500,  180.000 ),
    (  15, "CE2", "C", "sp2a",   13,   11,   10,    1.310,  109.000,    0.000 ),
    (  16, "CZ2", "C", "sp2a",   15,   13,   11,    1.400,  128.000,  180.000 ),
    (  17, "HZ2", "H", "",      16,   15,   13,    1.090,  120.000,    0.000 ),
    (  18, "CH2", "C", "sp2a",   16,   15,   13,    1.390,  116.000,  180.000 ),
    (  19, "HH2", "H", "",      18,   16,   15,    1.090,  120.000,  180.000 ),
    (  20, "CZ3", "C", "sp2a",   18,   16,   15,    1.350,  121.000,    0.000 ),
    (  21, "HZ3", "H",  "",     20,   18,   16,    1.090,  120.000,  180.000 ),
    (  22, "CE3", "C", "sp2a",   20,   18,   16,    1.410,  122.000,    0.000 ),
    (  23, "HE3", "H", "",      22,   20,   18,    1.090,  120.000,  180.000 ),
    (  24, "CD2", "C", "sp2a",   22,   20,   18,    1.400,  117.000,    0.000 ),
    (  25, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  26, "O  ", "O", "sp2",   25,    5,    3,    1.229,  120.500,    0.000 ),
]

TYR_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB2", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "HB3", "H", "",       7,    5,    3,    1.090,  109.500,   60.000 ),
    (  10, "CG ", "C", "sp2a",    7,    5,    3,    1.510,  109.470,  180.000 ),
    (  11, "CD1", "C", "sp2a",   10,    7,    5,    1.400,  120.000,  180.000 ),
    (  12, "HD1", "H", "",      11,   10,    7,    1.090,  120.000,    0.000 ),
    (  13, "CE1", "C", "sp2a",   11,   10,    7,    1.400,  120.000,  180.000 ),
    (  14, "HE1", "H", "",      13,   11,   10,    1.090,  120.000,  180.000 ),
    (  15, "CZ ", "C", "sp2a",   13,   11,   10,    1.400,  120.000,    0.000 ),
    (  16, "OH ", "O", "sp3",   15,   13,   11,    1.360,  120.000,  180.000 ),
    (  17, "HH ", "H", "",      16,   15,   13,    0.960,  113.000,    0.000 ),
    (  18, "CE2", "C", "sp2a",   15,   13,   11,    1.400,  120.000,    0.000 ),
    (  19, "HE2", "H", "",      18,   15,   13,    1.090,  120.000,  180.000 ),
    (  20, "CD2", "C", "sp2a",   18,   15,   13,    1.400,  120.000,    0.000 ),
    (  21, "HD2", "H", "",      20,   18,   15,    1.090,  120.000,  180.000 ),
    (  22, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  23, "O  ", "O", "sp2",   22,    5,    3,    1.229,  120.500,    0.000 ),
]

LYS_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB2", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "HB3", "H", "",       7,    5,    3,    1.090,  109.500,   60.000 ),
    (  10, "CG ", "C", "sp3",    7,    5,    3,    1.525,  109.470,  180.000 ),
    (  11, "HG2", "H", "",      10,    7,    5,    1.090,  109.500,  300.000 ),
    (  12, "HG3", "H", "",      10,    7,    5,    1.090,  109.500,   60.000 ),
    (  13, "CD ", "C", "sp3",   10,    7,    5,    1.525,  109.470,  180.000 ),
    (  14, "HD2", "H", "",      13,   10,    7,    1.090,  109.500,  300.000 ),
    (  15, "HD3", "H", "",      13,   10,    7,    1.090,  109.500,   60.000 ),
    (  16, "CE ", "C", "sp3",   13,   10,    7,    1.525,  109.470,  180.000 ),
    (  17, "HE2", "H", "",      16,   13,   10,    1.090,  109.500,  300.000 ),
    (  18, "HE3", "H", "",      16,   13,   10,    1.090,  109.500,   60.000 ),
    (  19, "NZ ", "N", "sp3",   16,   13,   10,    1.470,  109.470,  180.000 ),
    (  20, "HZ1", "H",  "",     19,   16,   13,    1.010,  109.470,   60.000 ),
    (  21, "HZ2", "H",  "",     19,   16,   13,    1.010,  109.470,  180.000 ),
    (  22, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  23, "O  ", "O", "sp2",   22,    5,    3,    1.229,  120.500,    0.000 ),
]

ARG_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB2", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "HB3", "H", "",       7,    5,    3,    1.090,  109.500,   60.000 ),
    (  10, "CG ", "C", "sp3",    7,    5,    3,    1.525,  109.470,  180.000 ),
    (  11, "HG2", "H", "",      10,    7,    5,    1.090,  109.500,  300.000 ),
    (  12, "HG3", "H", "",      10,    7,    5,    1.090,  109.500,   60.000 ),
    (  13, "CD ", "C", "sp3",   10,    7,    5,    1.525,  109.470,  180.000 ),
    (  14, "HD2", "H", "",      13,   10,    7,    1.090,  109.500,  300.000 ),
    (  15, "HD3", "H", "",      13,   10,    7,    1.090,  109.500,   60.000 ),
    (  16, "NE ", "N", "sp3",   13,   10,    7,    1.480,  111.000,  180.000 ),
    (  17, "HE ", "H", "",      16,   13,   10,    1.010,  118.500,    0.000 ),
    (  18, "CZ ", "C", "sp2",   16,   13,   10,    1.330,  123.000,  180.000 ),
    (  19, "NH1", "N", "sp3",   18,   16,   13,    1.330,  122.000,    0.000 ),
    (  20, "HH1", "H",  "",     19,   18,   16,    1.010,  119.800,    0.000 ),
    (  21, "HH1", "H",  "",     19,   18,   16,    1.010,  119.800,  180.000 ),
    (  22, "NH2", "N",  "sp2",  18,   16,   13,    1.330,  118.000,  180.000 ),
    (  23, "HH2", "H",  "",     22,   18,   16,    1.010,  119.800,    0.000 ),
    (  24, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  25, "O  ", "O", "sp2",   24,    5,    3,    1.229,  120.500,    0.000 ),
]

HIS_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB2", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "HB3", "H", "",       7,    5,    3,    1.090,  109.500,   60.000 ),
    (  10, "CG ", "C", "sp2",    7,    5,    3,    1.510,  115.000,  180.000 ),
    (  11, "ND1", "N", "sp3",   10,    7,    5,    1.390,  122.000,  180.000 ),
    (  12, "HD1", "H", "",      11,   10,    7,    1.010,  126.000,    0.000 ),
    (  13, "CE1", "C", "sp2",   11,   10,    7,    1.320,  108.000,  180.000 ),
    (  14, "HE1", "H", "",      13,   11,   10,    1.090,  120.000,  180.000 ),
    (  15, "NE2", "N", "sp2s",   13,   11,   10,    1.310,  109.000,    0.000 ),
    (  16, "CD2", "C", "sp2s",   15,   13,   11,    1.360,  110.000,    0.000 ),
    (  17, "HD2", "H", "",      16,   15,   13,    1.090,  120.000,  180.000 ),
    (  18, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  19, "O  ", "O", "sp2",   18,    5,    3,    1.229,  120.500,    0.000 ),
]

ASP_ZMATRIX = [
    (   0, "DUM",  "",  "",     -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM",  "",  "",      0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM",  "",  "",      1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ",  "N", "sp2(graphitic)",   2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ",  "H", "",      3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ",  "C", "sp3",   3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ",  "H", "",      5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ",  "C", "sp3",   5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB2",  "H", "",      7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "HB3",  "H", "",      7,    5,    3,    1.090,  109.500,   60.000 ),
    (  10, "CG ",  "C", "sp2",   7,    5,    3,    1.527,  109.470,  180.000 ),
    (  11, "OD1",  "O", "sp2",  10,    7,    5,    1.260,  117.200,   90.000 ),
    (  12, "OD2",  "O", "sp3",  10,    7,    5,    1.260,  117.200,  270.000 ),
    (  13, "HD2",  "H", "",     12,   10,    7,    0.960,  109.500,  180.000 ),
    (  14, "C  ",  "C", "sp2",   5,    3,    2,    1.522,  111.100,  180.000 ),
    (  15, "O  ",  "O", "sp2",  14,    5,    3,    1.229,  120.500,    0.000 ),
]

ASN_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB2", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "HB3", "H", "",       7,    5,    3,    1.090,  109.500,   60.000 ),
    (  10, "CG ", "C", "sp2",    7,    5,    3,    1.522,  111.100,  180.000 ),
    (  11, "OD1", "O", "sp2",   10,    7,    5,    1.229,  120.500,    0.000 ),
    (  12, "ND2", "N", "sp3",   10,    7,    5,    1.335,  116.600,  180.000 ),
    (  13, "HD2", "H", "",      12,   10,    7,    1.010,  119.800,  180.000 ),
    (  14, "HD2", "H", "",      12,   10,    7,    1.010,  119.800,    0.000 ),
    (  15, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  16, "O  ", "O", "sp2",   15,    5,    3,    1.229,  120.500,    0.000 ),
]

GLN_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "N", "sp2(graphitic)",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    1.010,  119.800,    0.000 ),
    (   5, "CA ", "C", "sp3",    3,    2,    1,    1.449,  121.900,  180.000 ),
    (   6, "HA ", "H", "",       5,    3,    2,    1.090,  109.500,  300.000 ),
    (   7, "CB ", "C", "sp3",    5,    3,    2,    1.525,  111.100,   60.000 ),
    (   8, "HB2", "H", "",       7,    5,    3,    1.090,  109.500,  300.000 ),
    (   9, "HB3", "H", "",       7,    5,    3,    1.090,  109.500,   60.000 ),
    (  10, "CG ", "C", "sp3",    7,    5,    3,    1.525,  109.470,  180.000 ),
    (  11, "HG2", "H", "",      10,    7,    5,    1.090,  109.500,  300.000 ),
    (  12, "HG3", "H", "",      10,    7,    5,    1.090,  109.500,   60.000 ),
    (  13, "CD ", "C", "sp2",   10,    7,    5,    1.522,  111.100,  180.000 ),
    (  14, "OE1", "O", "sp2",   13,   10,    7,    1.229,  120.500,    0.000 ),
    (  15, "NE2", "N", "sp3",   13,   10,    7,    1.335,  116.600,  180.000 ),
    (  16, "HE2", "H", "",      15,   13,   10,    1.010,  119.800,  180.000 ),
    (  17, "HE2", "H", "",      15,   13,   10,    1.010,  119.800,    0.000 ),
    (  18, "C  ", "C", "sp2",    5,    3,    2,    1.522,  111.100,  180.000 ),
    (  19, "O  ", "O", "sp2",   18,    5,    3,    1.229,  120.500,    0.000 ),
]

NTERM_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "N  ", "H", "",       2,    1,    0,    1.335,  116.600,  180.000 ),
]

# Note: this is just a fake "N  " label for the _buildResidue
# to make it thinking that we are starting a new amino acid.
CTERM_ZMATRIX = [
    (   0, "DUM", "",  "",      -1,   -2,   -3,    0.000,    0.000,    0.000 ),
    (   1, "DUM", "",  "",       0,   -1,   -2,    1.449,    0.000,    0.000 ),
    (   2, "DUM", "",  "",       1,    0,   -1,    1.522,  111.100,    0.000 ),
    (   3, "OXT", "O", "sp3",    2,    1,    0,    1.335,  116.600,  180.000 ),
    (   4, "H  ", "H", "",       3,    2,    1,    0.960,  109.500,    0.000 ),
]

# all amino acids, z-matrices and their sizes
AMINO_ACIDS = [
    ( "Alanine",       "ALA", "A", ALA_ZMATRIX, 13 ),
    ( "Arginine",      "ARG", "R", ARG_ZMATRIX, 26 ),
    ( "Asparagine",    "ARG", "N", ASN_ZMATRIX, 17 ),
    ( "Aspartic Acid", "ARG", "D", ASP_ZMATRIX, 16 ),
    ( "Cysteine",      "CYS", "C", CYS_ZMATRIX, 14 ),
    ( "Glutamic Acid", "GLU", "E", GLU_ZMATRIX, 19 ),
    ( "Glutamine",     "GLN", "Q", GLN_ZMATRIX, 20 ),
    ( "Glycine",       "GLY", "G", GLY_ZMATRIX, 10 ),
    ( "Histidine",     "HIS", "H", HIS_ZMATRIX, 20 ),
    ( "Isoleucine",    "ILE", "I", ILE_ZMATRIX, 22 ),
    ( "Leucine",       "LEU", "L", LEU_ZMATRIX, 22 ),
    ( "Lysine",        "LYS", "K", LYS_ZMATRIX, 24 ),
    ( "Methionine",    "MET", "M", MET_ZMATRIX, 20 ),
    ( "Phenylalanine", "PHE", "F", PHE_ZMATRIX, 23 ),
    ( "Proline",       "PRO", "P", PRO_ZMATRIX, 17 ),
    ( "Serine",        "SER", "S", SER_ZMATRIX, 14 ),
    ( "Threonine",     "THR", "T", THR_ZMATRIX, 17 ),
    ( "Tryptophan",    "TRP", "W", TRP_ZMATRIX, 27 ),
    ( "Tyrosine",      "TYR", "Y", TYR_ZMATRIX, 24 ),
    ( "Valine",        "VAL", "V", VAL_ZMATRIX, 19 )
]

# degrees to radians conversion
DEG2RAD = (pi/180.0)

def enablePeptideGenerator(enable):
    """
    This function enables/disables the Peptide Generator command by hiding or
    showing it in the Command Manager toolbar and menu.
    The enabling/disabling is done by the user via the "secret" NE1 debugging
    menu.

    To display the secret debugging menu, hold down Shift+Ctrl+Alt keys
    (or Shift+Cmd+Alt on Mac) and right click over the graphics area.
    Select "debug prefs submenu > Peptide Generator" and
    set the value to True. The "Peptide" option will then appear on the
    "Build" Command Manager toolbar/menu.

    @param enable: If true, the Peptide Generator is enabled. Specifically, it
                   will be added to the "Build" Command Manager toolbar and
                   menu.
    @type  enable: bool
    """
    win = env.mainwindow()
    win.buildProteinAction.setVisible(enable)

def get_unit_length(phi, psi):
    """
    Calculate a length of single amino acid in particular 
    secondary conformation.
    """
    # All unit length values obtained via measurements by me.
    # Fixes bug 2959. --Mark 2008-12-23
    if phi == -57.0 and psi == -47.0:
        unit_length = 1.5 # Alpha helix
    elif phi == -135.0 and psi == 135.0:
        unit_length = 3.4 # Beta strand
    elif phi == -55.0 and psi == -70.0:
        unit_length = 1.05 # Pi helix
    elif phi == -49.0 and psi == -26.0:
        unit_length = 1.95 # 3_10 helix
    elif phi == -75.0 and psi == 150.0:
        unit_length = 3.14 # Polyproline-II helix
    elif phi == -180.0 and psi == 180.0:
        unit_length = 3.6 # Fully extended
    else:
        # User chose "Custom" conformation option in the Insert Peptide PM 
        # which lets the user set any phi-psi angle values.
        # We need a formula to estimate the proper unit_length given the
        # conformational angles phi and psi. It would also be a good idea
        # to check these values to confirm they are an allowed combination.
        # For more info, do a google search for "Ramachandran plot".
        # For now, I'm setting the unit length to 1.5. --Mark 2008-12-23
        unit_length = 1.5
        msg = "\nUser selected custom conformational angles: "\
            "phi=%.2f, psi=%.2f.\nSetting unit_length=%.2f\n" % \
            (phi, psi, unit_length)
        print_compact_stack(msg)
        
    return unit_length
    
class PeptideGenerator:
    prev_coords = zeros([3,3], Float)

    peptide_mol = None
    length = 0
    prev_psi = 0
    
    # Based on analogous Nanotube Builder method.
    def _orient(self, chunk, pt1, pt2):
        """
        Orients the Peptide I{chunk} based on two points. I{pt1} is
        the first endpoint (origin) of the Peptide. The vector I{pt1}, I{pt2}
        defines the direction and central axis of the Peptide.
        
        piotr 080801: I copied this method from Nanotube Builder.
        
        @param pt1: The starting endpoint (origin) of the Peptide.
        @type  pt1: L{V}
        
        @param pt2: The second point of a vector defining the direction
                    and central axis of the Peptide.
        @type  pt2: L{V}
        """
        
        a = V(0.0, 0.0, -1.0)
        # <a> is the unit vector pointing down the center axis of the default
        # structure which is aligned along the Z axis.
        bLine = pt2 - pt1
        bLength = vlen(bLine)
        if bLength == 0:
            return
        b = bLine / bLength
        # <b> is the unit vector parallel to the line (i.e. pt1, pt2).
        axis = cross(a, b)
        # <axis> is the axis of rotation.
        theta = angleBetween(a, b)
        # <theta> is the angle (in degress) to rotate about <axis>.
        scalar = bLength * 0.5
        rawOffset = b * scalar
        
        if theta == 0.0 or theta == 180.0:
            axis = V(0, 1, 0)
            # print "Now cross(a,b) =", axis
            
        rot =  (pi / 180.0) * theta  # Convert to radians
        qrot = Q(axis, rot) # Quat for rotation delta.
        
        # Move and rotate the Peptide into final orientation.
        
        chunk.move(-chunk.center)

        chunk.rot(qrot)
        
        # Bruce suggested I add this. It works here, but not if its 
        # before move() and rot() above. Mark 2008-04-11
        chunk.full_inval_and_update()
        return
        
    def get_number_of_res(self, pos1, pos2, phi, psi):
        """
        Calculate a number of residues necessary to fill 
        the pos1-pos2 vector.
        
        @param pos1, pos2: vector points
        @type pos1, pos2: V
        
        @param phi, psi: peptide chain angles
        @type phi, psi: float
        """
        return 1 + int(vlen(pos2 - pos1) / get_unit_length(phi, psi))
    
    def make_aligned(self, 
                     assy, 
                     name, 
                     aa_idx, 
                     phi, psi, 
                     pos1, pos2, 
                     secondary = SS_COIL, 
                     fake_chain = False, 
                     length = None):
        """
        Build and return a chunk that is a homo-peptide aligned to 
        a pos2-pos1 vector.
        
        @param aa_idx: amino acid type (index in AMINO_ACIDS list)
        @type aa_idx: int
        
        @param name: chunk name
        @type name: string
        
        @param phi, psi: peptide bond angles 
        @type phi, psi: float
        
        @param pos1, pos2: desired peptide positions (beginning and end)
        @type pos1, pos2: V
        
        @param secondary: secondary structure class, used for visual representation
        The actual peptide chain conformation is based on phi / psi angles.
        @type secondary: int
        
        @param fake_chain: if True, create only C-alpha atoms. used for drawing
        peptide trace image during interactive peptide placement (used by
        PeptideLine_GraphicsMode.py)
        @type fake_chain: boolean
        
        @param length: optional peptide length (number of amino acids), if 
        not specified, pos1 and pos2 are used to figure out the length
        @type length: int
        
        @return: A homo-polypeptide chain.
        @rtype:  L{Chunk}
        """

        if not length:
            self.length = self.get_number_of_res(pos1, pos2, phi, psi)
            if self.length == 0:
                return None
        else:
            self.length = length
            
        # Create a molecule
        mol = Chunk(assy, name)
            
        if not fake_chain:
            mol.protein = Protein()
            mol.protein.set_chain_id('A')
            
        # Generate dummy atoms positions
        self.prev_coords[0][0] = pos1[0] - 1.0
        self.prev_coords[0][1] = pos1[1] - 1.0
        self.prev_coords[0][2] = pos1[2] 

        self.prev_coords[1][0] = pos1[0] - 1.0
        self.prev_coords[1][1] = pos1[1] 
        self.prev_coords[1][2] = pos1[2] 

        self.prev_coords[2][0] = pos1[0]
        self.prev_coords[2][1] = pos1[1]
        self.prev_coords[2][2] = pos1[2]

        name, short_name, symbol, zmatrix, size = AMINO_ACIDS[aa_idx]
        
        # Add a N-terminal hydrogen
        self.nterm_hydrogen = None

        # Initially, the Peptide Builder was creating peptide structures
        # saturated at both ends, i.e. with N-terminal hydrogen and C-terminal
        # OH group present. Currently, this code is commented out to allow 
        # connecting multiple peptide structure be creating bonds between
        # the C- and N- terminal ends of two individual structures.
        """
        if not fake_chain:
            atom = Atom("H", pos1, mol)
            atom._is_aromatic = False
            atom._is_single = False
            self.nterm_hydrogen = atom
            mol.protein.add_pdb_atom(atom, "H", 1, name)
            atom.pdb_info = {}
            atom.pdb_info['atom_name'] = "H"
            atom.pdb_info['residue_name'] = short_name
            atom.pdb_info['residue_id'] = "  1 "
            atom.pdb_info['standard_atom'] = True
        """
    
        self.init_ca = None
        
        # Generate the peptide chain.
        for idx in range(int(self.length)):
            self._buildResidue(mol, zmatrix, size, idx+1, phi, psi, secondary, None, short_name, fake_chain=fake_chain)

        # See the comment above.
        """
        # Add a C-terminal OH group
        self._buildResidue(mol, CTERM_ZMATRIX, 5, int(self.length), 0.0, 0.0, secondary, None, short_name, fake_chain=fake_chain)        
        """
        
        # Compute bonds (slow!)
        # This should be replaced by a proper bond assignment.
        
        if not fake_chain:
            inferBonds(mol)

        # Assign proper bond orders.
        i = 1
        for atom in mol.atoms.itervalues():
            if atom.bonds:
                for bond in atom.bonds:
                    if bond.atom1.getAtomTypeName()=="sp2" and \
                       bond.atom2.getAtomTypeName()=="sp2":
                        if (bond.atom1._is_aromatic and
                            bond.atom2._is_aromatic):
                            bond.set_v6(V_AROMATIC)
                        elif ((bond.atom1._is_aromatic == False and
                               bond.atom1._is_aromatic == False) and
                               not (bond.atom1._is_single and
                                    bond.atom2._is_single)):
                            bond.set_v6(V_DOUBLE)
            i += 1
                            
        # Remove temporary attributes.
        for atom in mol.atoms.itervalues():
            del atom._is_aromatic
            del atom._is_single

        # Axis of first selected chunk
        ax = V(0.,0.,1.) 
        mol.rot(Q(mol.getaxis(),ax))
        
        self._orient(mol, pos2, pos1)
        
        if self.init_ca:
            mol.move(pos1 - self.init_ca.posn())
        
        mol_dummy = None
        
        return mol          

    def _buildResidue(self, mol, zmatrix, n_atoms, idx, phi, psi, secondary, init_pos, residue_name, fake_chain=False):
        """
        Builds cartesian coordinates for an amino acid from the internal
        coordinates table.

        @param mol: a chunk to which the amino acid will be added.
        @type mol: Chunk
        
        @param zmatrix: is an internal coordinates array corresponding to a 
        given amino acid.
        @type zmatrix: list
        
        @param n_atoms: size of z-matrix (a number of atoms to be build + 3 
        dummy atoms)
        @type n_atoms: int

        @param idx: is a residue index (1..length).
        @type idx: integer
        
        @param phi, psi: peptide bond phi and psi angles
        @type phi, psi: float
        
        @param init_pos: optional postions of previous CA, C and O atoms.
        @type init_pos: V
        
        @param symbol: current amino acid symbol (used to derermine proline case)
        @type symbol: string
        
        """

        # note: currently, it doesn't rebuild bonds, so inferBonds has to be 
        # called after this method. Unfortunately, the proper bond order can 
        # not be correctly recognized this way. Therefore, temporary atom flags 
        # _is_aromatic and _is_single are used.
        
        #this code was re-factored by EricM and internal-to-cartesian
        # conversion method was moved to geometry.InternalCoordinatesToCartesian

        if mol is None:
            return

        if not init_pos: # assign three previous atom positions
            coords = self.prev_coords
        else: 
            # if no prev_coords are given, compute the first three atom positions
            coords = zeros([3,3], Float)
            num, name, atom_name, atom_type, \
               atom_c, atom_b, atom_a, r, a, t = zmatrix[1]
            coords[0][0] = 0.0;
            coords[0][1] = 0.0;
            coords[0][2] = 0.0;
            coords[1][0] = r;
            coords[1][1] = 0.0;
            coords[1][2] = 0.0;
            ccos = cos(DEG2RAD*a)
            num, name, atom_name, atom_type, \
               atom_c, atom_b, atom_a, r, a, t = zmatrix[2]
            if atom_c == 1:
                coords[2][0] = coords[0][0] + r*ccos
            else:
                coords[2][0] = coords[0][0] - r*ccos
            coords[2][1] = r * sin(DEG2RAD*a)
            coords[2][2] = 0.0
            for i in range (0, 3):
                self.prev_coords[i][0] = coords[i][0] + init_pos[0]
                self.prev_coords[i][1] = coords[i][1] + init_pos[1]
                self.prev_coords[i][2] = coords[i][2] + init_pos[2]

        translator = InternalCoordinatesToCartesian(n_atoms, coords)

        for n in range (3, n_atoms):
            # Generate all coordinates using three previous atoms
            # as a frame of reference,
            num, name, atom_name, atom_type, \
               atom_c, atom_b, atom_a, r, a, t = zmatrix[n]

            # Apply the peptide bond conformation
            if residue_name != "PRO":
                if name == "N  " and not init_pos:
                    t = self.prev_psi + 0.0
                if name == "O  ":
                    t = psi + 180.0
                if name == "HA " or name == "HA2":
                    t = 120.0 + phi
                if name == "CB " or name == "HA3":
                    t = 240.0 + phi
                if name == "C  ":
                    t = phi
            else:
                # proline
                if name == "N  " and not init_pos:
                    t = self.prev_psi + 0.0
                if name == "O  ":
                    t = psi + 180.0
                if name == "CA ":
                    t = phi - 120.0
                if name == "CD ":
                    t = phi + 60.0

            translator.addInternal(n+1, atom_c+1, atom_b+1, atom_a+1, r, a, t)
            xyz = translator.getCartesian(n+1)

            if self.nterm_hydrogen:
                # This is a hack for the first hydrogen atom
                # to make sure the bond length is correct.
                self.nterm_hydrogen.setposn(
                    self.nterm_hydrogen.posn() + \
                    0.325 * norm(xyz))
                self.nterm_hydrogen = None

            # Store previous coordinates for the next building step
            if not init_pos:
                if name=="N  ":
                    self.prev_coords[0][0] = xyz[0]
                    self.prev_coords[0][1] = xyz[1]
                    self.prev_coords[0][2] = xyz[2]
                if name=="CA ":
                    self.prev_coords[1][0] = xyz[0]
                    self.prev_coords[1][1] = xyz[1]
                    self.prev_coords[1][2] = xyz[2]
                if name=="C  ":
                    self.prev_coords[2][0] = xyz[0]
                    self.prev_coords[2][1] = xyz[1]
                    self.prev_coords[2][2] = xyz[2]
                    
            # Add a new atom to the molecule
            if not fake_chain or \
               name == "CA ":
                atom = Atom(
                    atom_name,
                    xyz,
                    mol)

                if not self.init_ca and \
                   name == "CA ":
                    self.init_ca = atom

                if mol.protein:
                    aa = mol.protein.add_pdb_atom(atom, 
                                             name.replace(' ',''), 
                                             idx, 
                                             AA_3_TO_1[residue_name])
                    atom.pdb_info = {}
                    atom.pdb_info['atom_name'] = name.replace(' ','')
                    atom.pdb_info['residue_name'] = residue_name
                    residue_id = "%3d " % idx
                    atom.pdb_info['residue_id'] = residue_id
                    atom.pdb_info['standard_atom'] = True
                    atom.pdb_info['chain_id'] = True
                    if aa:
                        aa.set_secondary_structure(secondary)

                # Create temporary attributes for proper bond assignment.
                atom._is_aromatic = False
                atom._is_single = False

                if atom_type == "sp2a":
                    atom_type = "sp2"
                    atom._is_aromatic = True

                if atom_type == "sp2s":
                    atom_type = "sp2"
                    atom._is_single = True

                atom.set_atomtype_but_dont_revise_singlets(atom_type)
                
                ### debug - output in PDB format	
                ### print "ATOM  %5d  %-3s %3s %c%4d    %8.3f%8.3f%8.3f" % ( n, name, "ALA", ' ', res_num, xyz[0], xyz[1], xyz[2])	

        self.prev_psi = psi # Remember previous psi angle.

        return
    
# end