! generate_water.cns ! soaks a protein structure in a layer of water ! can be applied iteratively (using dyncount > 1) ! ! ************************************ ! * Authors and copyright: * ! * Michael Nilges, Jens Linge, EMBL * ! * No warranty implied or expressed * ! * All rights reserved * ! ************************************ eval ($boxlength = 18.856) ! length of Brooks' water box eval ($thickness = 8) ! maxi. initial water-protein distance (heavy atoms) eval ($pw_dist = 4.0) ! mini. initial water-protein distance (heavy atoms) eval ($water_diam = 2.4) ! diameter of water molecule eval ($dyncount = 1) ! iteration number (usually 1) eval ($water = "WAT" + encode($dyncount)) !-------------------------------------------------- ! read in the same box of water several times, and move it around ! so as to cover all the space around the site of interest. ! take into account box offset show max (x) ((segid PRO*)) evaluate ($xmax = $result) show min (x) ((segid PRO*)) evaluate ($xmin = $result) show max (y) ((segid PRO*)) evaluate ($ymax = $result) show min (y) ((segid PRO*)) evaluate ($ymin = $result) show max (z) ((segid PRO*)) evaluate ($zmax = $result) show min (z) ((segid PRO*)) evaluate ($zmin = $result) ! loop over several iterations of water filling and dynamics !-------------------------------------------------- ! read in the same box of water several times, and move it around ! so as to cover all the space around the site of interest. ! take into account box offset ! determine how many boxes are necessary in each dimension eval ($xbox = int( ($xmax - $xmin + 2 * ($thickness + $water_diam)) / $boxlength + 0.5)) eval ($ybox = int( ($ymax - $ymin + 2 * ($thickness + $water_diam)) / $boxlength + 0.5)) eval ($zbox = int( ($zmax - $zmin + 2 * ($thickness + $water_diam)) / $boxlength + 0.5)) eval ($xmtran = $xmax + $thickness - $boxlength/2 + $water_diam) eval ($ymtran = $ymax + $thickness - $boxlength/2 + $water_diam) eval ($zmtran = $zmax + $thickness - $boxlength/2 + $water_diam) set echo off message off end eval ($xcount=0) eval ($xtrans = $xmin - $thickness - $water_diam - $boxlength/2 ) while ($xtrans < $xmtran) loop wat1 eval ($xcount=$xcount+1) eval ($xtrans = $xtrans + $boxlength) eval ($ycount=0) eval ($ytrans = $ymin - $thickness - $water_diam - $boxlength/2 ) while ($ytrans < $ymtran) loop wat2 eval ($ycount=$ycount+1) eval ($ytrans = $ytrans + $boxlength) eval ($zcount=0) eval ($ztrans = $zmin - $thickness - $water_diam - $boxlength/2 ) while ($ztrans < $zmtran) loop wat3 eval ($zcount=$zcount+1) eval ($ztrans = $ztrans + $boxlength) segment name=" " chain coordinates @@boxtyp20.pdb end end coor @@boxtyp20.pdb do (segid=W000) (segid " ") coor sele=(segid W000) translate vector = ($xtrans $ytrans $ztrans) end ident (store1) (segid W000 and name oh2) ident (store2) (store1 and ((segid PRO*) and not hydro) around $pw_dist) ident (store3) (store1 and (segid wat# and not hydro) around $water_diam) ident (store4) (store1 and not ((segid PRO*) and not hydro) around $thickness) delete sele= (byres (store2 or store3 or store4)) end ! give waters unique segid name eval ($segid= "W" + encode($xcount) + encode($ycount) + encode($zcount)) do (segid = $segid) (segid W000) end loop wat3 end loop wat2 end loop wat1 ! now, give waters a unique resid so that we get the segid to play around with ident (store1) (all) show min (store1) (segid w*) do (store1 = store1 - $result + 1) (segid w*) do (resid = encode(int(store1/3 -0.1) +1)) (segid w* and not segid wat#) do (segid = $water) (segid w* and not segid wat#) ! shave off any waters that left delete sele= (byres (name oh2 and not ((segid PRO*) and not hydro) around $thickness)) end {* write out initial coordinates *} evaluate ($filename=$fileroot+ encode($count)+ "wini.pdb") ! write coordinates sele= (all) output =$filename end set echo on message on end