remarks   PHI/PSI plot
remarks   Kleywegt/Jones contours (Structure 4:1395-1400, 1997)
remarks	  Ethan A Merritt - March 1998

structure @ generate.psf end

parameter @ TOPPAR:parhcsdx.pro end

eval ($CHAINID_IS_SEGID = 1)

 coordinates           @ current.protein

delete select=(not known) end


set display=rama.gnu end             ! write input file for gnuplot

!
! The Symbol font is required to get the axis labels and the
! degree signs.  Unfortunately this means that the title will be
! printed in Greek.  To fix this you can eventually edit the
! PostScript file to change the title font to something more reasonable.
! (That basically means inserting a single line)
!

display set term postscript portrait mono solid "Symbol" 10
display set output 'rama.ps'

display set title 'Ramachandran plot $DATE '
display set xlabel "F"
display set xtic   -180., 60., 180.
display set xrange [-180:180]
display set ylabel "Y"
display set ytic   -180., 60., 180.
display set yrange [-180:180]
display set format xy "%4.0f°"
display set border
display set size 0.7, 1.0
display set data style lines
 
set echo=false end

!
! You must tailor this to your structure
! The selection must run from the 2nd to the 2nd-last residue of each chain
!
vector ident (store4)
	((segid D or segid E or segid F or segid G or segid H) and
	 ( name CA and resid 2:102 ))


for $1 in id ( recall 4 ) loop main

    vector show element (segid) (id $1)
    evaluate ($SEG = $result)
    vector show element (resid) (id $1)
    evaluate ($3 = $result)
    evaluate ($2 = DECODE($3) - 1)
    evaluate ($4 = DECODE($3) + 1)
    pick dihedral
	(SEGID $SEG and resid $2 and name C)
	(SEGID $SEG and resid $3 and name N)
	(SEGID $SEG and resid $3 and name CA)
	(SEGID $SEG and resid $3 and name C)
	geometry
    evaluate ($5 = INT($RESULT))
    pick dihedral
	(SEGID $SEG and resid $3 and name N)
	(SEGID $SEG and resid $3 and name CA)
	(SEGID $SEG and resid $3 and name C)
	(SEGID $SEG and resid $4 and name N)
	geometry
    evaluate ($6 = INT($RESULT))
!
! Uncomment this one if you want the residue number plotted,
!   display set label ' $3 ' at $5, $6 center
!
! Otherwist Plot a + for glycine; a dot for all other residue types
!
    vector show element (resname) (id $1)
    evaluate ($resname = $result)
    if ($resname = "GLY") then
   	display set label '+' at $5, $6 center
    else
   	display set label '·' at $5, $6 center
    end if


end loop main

display set nokey
display set nosurface
display set contour base
display set parametric
display set cntrparam levels discrete 20., 225.
display set cntrparam bspline
display set view 0, 0

display splot 'kleywegt_matrix.gnu' with line 2


set echo=true end
stop