charmm_prm_to_rosetta_mm.ipynb

{
 "worksheets": [
  {
   "cells": [
    {
     "metadata": {},
     "cell_type": "code",
     "input": "\"\"\"\nCreated on Tue Sep 30 11:04:24 2014\n\nA script to convert Charmm parameter files into mm_atom_properties.txt and\nmm_torsion_params.txt used by Rosetta.\n\nUses ParamEd (https://github.com/swails/ParmEd)\n\nIdeally this would merge/update the existing Rosetta parameters. \nFor now the resulting parameters only works (have only been tested) in slim mode, meaning\nonly proteins (and the R1A spin label)\n@author: ajasja\n\"\"\";",
     "prompt_number": 7,
     "outputs": [],
     "language": "python",
     "trusted": true,
     "collapsed": false
    },
    {
     "metadata": {},
     "cell_type": "code",
     "input": "import glob\nimport chemistry\nimport chemistry.charmm.parameters as pr\nimport zlib\nimport datetime\nprint \"ParamEd version\", chemistry.__version__\n\ndef get_crc(name):\n     return hex(zlib.crc32(open(name).read()) & 0xffffffff)\n",
     "prompt_number": 8,
     "outputs": [
      {
       "output_type": "stream",
       "text": "ParamEd version 1.1.1\n",
       "stream": "stdout"
      }
     ],
     "language": "python",
     "trusted": true,
     "collapsed": false
    },
    {
     "metadata": {},
     "cell_type": "code",
     "input": "#chose files to convert\nrtf_files = glob.glob('../FF/*.rtf')\nprm_files = glob.glob('../FF/*.prm')\nprm_crc =[get_crc(name) for name in prm_files]\nprm_files",
     "prompt_number": 9,
     "outputs": [
      {
       "text": "['../FF\\\\par_all27_prot_lipid.prm', '../FF\\\\par_all27_R1A.prm']",
       "output_type": "pyout",
       "metadata": {},
       "prompt_number": 9
      }
     ],
     "language": "python",
     "trusted": true,
     "collapsed": false
    },
    {
     "metadata": {},
     "cell_type": "code",
     "input": "#read in the parameter set\nparam_set = pr.CharmmParameterSet()\nfor rtf_file in rtf_files:\n    param_set.read_topology_file(rtf_file)\nfor prm_file in prm_files:\n    param_set.read_parameter_file(prm_file)\n    \nparam_set.condense()",
     "prompt_number": 10,
     "outputs": [
      {
       "text": "<chemistry.charmm.parameters.CharmmParameterSet at 0x6e386b0>",
       "output_type": "pyout",
       "metadata": {},
       "prompt_number": 10
      }
     ],
     "language": "python",
     "trusted": true,
     "collapsed": false
    },
    {
     "metadata": {},
     "cell_type": "code",
     "input": "rosetta_specific_atom_types =\"\"\"\\\nON3\t-0.120000\t1.700000\t-0.120000\t1.700000\nP\t-0.585000\t2.150000\t-0.585000\t2.150000\nB2\t-0.200000\t1.980000\t-0.200000\t1.980000 # Copied from NR1 with modification to the radius  \nBR\t-0.320000\t2.220000\t-0.320000\t2.220000 # radius/wdepth: taken from MM3 forcefield\nI\t-0.424000\t2.360000\t-0.424000\t2.360000 # radius/wdepth: taken from MM3 forcefield\nCO\t-0.250000\t1.568000\t-0.250000\t1.568000 # cobalt, modified from zinc\nCU\t-0.250000\t1.097800\t-0.250000\t1.097800 # copper, modified from zinc\nX\t 0.000000\t0.000000\t 0.000000\t0.000000\nVIRT\t 0.000000\t0.000000\t 0.000000\t0.000000\n\"\"\"\n\n#write atom properties\nwith open('mm_atom_properties.txt', 'w') as f:\n    #f.write('#Generated from: ' + str(prm_files)+'\\n')\n    #f.write('#CRC32 of files:' + str(prm_crc)+'\\n')\n    #f.write('#'+str(datetime.datetime.today())+'\\n')\n    #Hmm atom properties can not have a header with comments?\n    \n    f.write(\"NAME LJ_WDEPTH LJ_RADIUS LJ_3B_WDEPTH LJ_3B_RADIUS\\n\")\n    fmt_str = \"{:<5s} {: 1.6f} {: 1.6f} {: 1.6f} {: 1.6f}\\n\"\n    for name, atom in param_set.atom_types_str.iteritems():        \n        #print name, atom.epsilon, atom.rmin, atom.epsilon_14, atom.rmin_14    \n        if atom.epsilon is None:\n            #some virtual particles don't have LJ info.\n            out_str = fmt_str.format(name, 0, 0, 0, 0)        \n        else:        \n            out_str =  fmt_str.format(name, atom.epsilon, atom.rmin, atom.epsilon_14, atom.rmin_14)\n        f.write(out_str)  \n    \n    f.write(rosetta_specific_atom_types)\n\n \n\n        \n        \n!head mm_atom_properties.txt    \n!tail mm_atom_properties.txt\nimport shutil\nshutil.copy('mm_atom_properties.txt', \n            '\\\\\\\\mint\\\\ajasja\\\\PyRosetta.namespace-r21\\\\database\\\\chemical\\\\mm_atom_type_sets\\\\fa_standard')",
     "prompt_number": 11,
     "outputs": [
      {
       "output_type": "stream",
       "text": "NAME LJ_WDEPTH LJ_RADIUS LJ_3B_WDEPTH LJ_3B_RADIUS\nH     -0.046000  0.224500 -0.046000  0.224500\nHC    -0.046000  0.224500 -0.046000  0.224500\nHA    -0.022000  1.320000 -0.022000  1.320000\nHT    -0.046000  0.224500 -0.046000  0.224500\nHP    -0.030000  1.358200 -0.030000  1.358200\nHB    -0.022000  1.320000 -0.022000  1.320000\nHR1   -0.046000  0.900000 -0.046000  0.900000\nHR2   -0.046000  0.700000 -0.046000  0.700000\nHR3   -0.007800  1.468000 -0.007800  1.468000\n",
       "stream": "stdout"
      },
      {
       "output_type": "stream",
       "text": "NN    -0.200000  1.850000 -0.200000  1.850000\nON3\t-0.120000\t1.700000\t-0.120000\t1.700000\nP\t-0.585000\t2.150000\t-0.585000\t2.150000\nB2\t-0.200000\t1.980000\t-0.200000\t1.980000 # Copied from NR1 with modification to the radius  \nBR\t-0.320000\t2.220000\t-0.320000\t2.220000 # radius/wdepth: taken from MM3 forcefield\nI\t-0.424000\t2.360000\t-0.424000\t2.360000 # radius/wdepth: taken from MM3 forcefield\nCO\t-0.250000\t1.568000\t-0.250000\t1.568000 # cobalt, modified from zinc\nCU\t-0.250000\t1.097800\t-0.250000\t1.097800 # copper, modified from zinc\nX\t 0.000000\t0.000000\t 0.000000\t0.000000\nVIRT\t 0.000000\t0.000000\t 0.000000\t0.000000\n",
       "stream": "stdout"
      }
     ],
     "language": "python",
     "trusted": true,
     "collapsed": false
    },
    {
     "metadata": {},
     "cell_type": "code",
     "input": "#write torsions\nwith open('mm_torsion_params.txt', 'w') as f:\n    #f.write('#Generated from: ' + str(prm_files)+'\\n')\n    #f.write('#CRC32 of files:' + str(prm_crc)+'\\n')\n    #f.write('#'+str(datetime.datetime.today())+'\\n')\n    f.write(\"#atom1 atom2 atom3 atom4 strength periodicity minimum\\n\")\n    fmt_str = \"{:<5s} {:<5s} {:<5s} {:<5s} {: >9.6f} {:1d} {: >7.3f}\\n\"\n    for names, diheds in param_set.dihedral_types.iteritems():            \n        for dihed in diheds:\n            #print names, dihed\n            f.write(fmt_str.format(names[0], names[1], names[2], names[3], \n                                 dihed.phi_k, dihed.per, dihed.phase) )\n            \n!head mm_torsion_params.txt           \n!tail mm_torsion_params.txt ",
     "prompt_number": 12,
     "outputs": [
      {
       "output_type": "stream",
       "text": "#atom1 atom2 atom3 atom4 strength periodicity minimum\nNH2   CT1   C     O      0.000000 1   0.000\nNH1   C     CT1   NH2    0.000000 1   0.000\nCT1   CT1   NH2   H      0.000000 1   0.000\nC     CT1   NH2   H      0.000000 1   0.000\nH     NH2   CT1   HB     0.110000 3   0.000\nCT2   CT1   NH2   H      0.110000 3   0.000\nCT3   CT1   NH2   H      0.110000 3   0.000\nX     FE    SS    X      0.000000 4   0.000\nX     CS    SS    X      0.000000 3   0.200\n",
       "stream": "stdout"
      },
      {
       "output_type": "stream",
       "text": "NH1   CT1   CT2   SM     0.237000 3   0.000\nC     CT1   CT2   SM     0.233000 1   0.000\nC     CT1   CT2   SM    -0.005000 2   0.000\nC     CT1   CT2   SM     0.149000 3   0.000\nHAL1  CTL1  CTL2  SM     0.064000 3   0.000\nHB    CT1   CT2   SM     0.064000 3   0.000\nCTL3  CTL1  CTL2  HAL2   0.120000 3   0.000\nCTL2  CEL1  CEL1  HEL1   1.000000 2 180.000\nCTL1  CEL1  CEL1  HEL1   1.000000 2 180.000\nCTL3  CTL1  CEL1  HEL1   0.120000 3   0.000\n",
       "stream": "stdout"
      }
     ],
     "language": "python",
     "trusted": true,
     "collapsed": false
    }
   ],
   "metadata": {}
  }
 ],
 "metadata": {
  "name": "",
  "signature": "sha256:b2568433ff1c8bbd2eeac105f1619c04c7f921cd06c26456691a02c076e2a43b",
  "gist_id": "0f62fc1afe23da134893"
 },
 "nbformat": 3
}