#
# Load settings
#

source 1CK7_settings.bcl

#
# Create non-standard residue library
#

createStdGroup NAME

copyStdResidue aminoAcids94/HID NAME/HD1
copyStdResidue aminoAcids94/HID NAME/HD2
copyStdResidue aminoAcids94/HID NAME/HD3
copyStdResidue aminoAcids94/CYM NAME/CY1
copyStdResidue metals/.ZN NAME/ZN1

#
# Create new atom types
#

# copy sulfure at the ruduced state (Cornell et al, 1995)
# SH amber atom type is copied rather then S, otherwise you will get an error using gaussian in further stages.
# this step is important for the assignment of the right valence and connectivity
copyAtomType parm99/SH NAME/SA 

# periods are replaced by spaces internally
setAtomType NAME/CY1/.SG. NAME/SA

#
# write standard library
#
writeLib NAME NAME.xml

#
# Open PDB file
#
readPdb NAME NAME_fixed.pdb

# Assign S-S bonds
assignDisulfideBonds

# Atom type
atomType

# Assign bonds, angles, and torsions
assignConnectivity
# Add hydrogens
addHs /NAME/

#
# Build Cluster
#
# Create molecule named CLR
createMolecule CLR

# Create HD1-1
createResidue HD1 in CLR
addToResidue /NAME/CLR/HD1 /NAME/1/HID-407 not bb

# Create HD1-2
createResidue HD2 in CLR
addToResidue /NAME/CLR/HD2 /NAME/1/HID-413 not bb

# Create HD1-3
createResidue HD3 in CLR
addToResidue /NAME/CLR/HD3 /NAME/1/HID-403 not bb

# Create CY1-4
createResidue CY1 in CLR  
addToResidue /NAME/CLR/CY1 /NAME/1/CYM-102 not bb

# Create ZN1-5 from 
createResidue ZN1 in CLR
addToResidue /NAME/CLR/ZN1 /NAME/2/ZN.-1000/ZN..

# Methyl terminating groups (Add these after all other residue have been added)
addFragment terminal/CH3 bd /NAME/CLR/HD1-1/.CB. ag /NAME/CLR/HD1-1/.CG. tr /NAME/CLR/HD1-1/.ND1 -57.5
addFragment terminal/CH3 bd /NAME/CLR/HD2-2/.CB. ag /NAME/CLR/HD2-2/.CG. tr /NAME/CLR/HD2-2/.ND1 153.3
addFragment terminal/CH3 bd /NAME/CLR/HD3-3/.CB. ag /NAME/CLR/HD3-3/.CG. tr /NAME/CLR/HD3-3/.ND1 117.7
addFragment terminal/CH3 bd /NAME/CLR/CY1-4/.CB. ag /NAME/CLR/CY1-4/.SG. tr /NAME/CLR/ZN1-5/ZN.. -169.4


# Create bonds with zinc1
createBond /NAME/CLR/ZN1/ZN.. /NAME/CLR/HD1/.NE2
createBond /NAME/CLR/ZN1/ZN.. /NAME/CLR/HD2/.NE2
createBond /NAME/CLR/ZN1/ZN.. /NAME/CLR/HD3/.NE2
createBond /NAME/CLR/ZN1/ZN.. /NAME/CLR/CY1/.SH.
# Atomtype
atomType /NAME/CLR
# Assign bonds and angles
assignParameters /NAME/CLR

# Add Bond and Angle Parameters
addBondAndAngleParameters /NAME/CLR NAME
# Write parameters file
writeParams NAME NAME_params.xml

# Write pdb/mol files (Good idea to view one of these files
# before running the Gaussian calc.)
writePdb /NAME/CLR NAME_sidechain.pdb
writeMol /NAME/CLR NAME_sidechain.mol
writeSdf /NAME/CLR NAME_sidechain.sdf

# Gaussian options
levelOfTheory B3LYP
basisSet 6-31G*
clusterCharge CLR 0
clusterSpin 1
gaussianNProc 4
gaussianMem 2GB
# Set Gaussian input name --> Optimize and get force constants
gaussianOptAndFC /NAME/CLR NAME_sidechain.com
# Exit MCPB
quit