# Copyright (c) 2002-08 Peter Guntert. All rights reserved.
## 7MACROS: ramaplot - CYANA macro
##
## Parameters: file=<filename>.ps             (default: ramaplot.ps)
##             nobackground
##             label
##
## Produces a graphics output file with the given name (a Graf file if the
## extension is ".grf", or a Postscript file otherwise) with a Ramachandran plot 
## of the selected structures. The background consists of three different blue
## tones indicating "most favored regions" (dark blue), "additional allowed
## regions" (medium blue) and "generously allowed regions" (light blue). It
## corresponds to the background found in the program Procheck (Laskowski et
## al., 1993).
## Optionally, the background can be omitted (nobackground) and the residues
## with backbone angles that lie outside of the allowed regions are labeled.

var echo c0 c1 c2 c3 phi psi x y xx yy i j k lastc lastx fillp str struct
var p1 p2 p3 p4 p5 p6 p7 p8 p9 p10 p11 p12 p13 p14 p15 p16 p17 p18 p19 p20
var p21 p22 p23 p24 p25 p26 p27 p28 p29 p30 p31 p32 p33 p34 p35 p36
var aas rp r rn s0 s1 s2 s3 n

syntax $macro:file=@f.ps=ramaplot.ps nobackground label

echo:=off

p1 := 2,3,3,3,3,3,3,3,3,3, 3,2,2,2,1,1,0,0,0,0 ,0,1,1,1,1,1,0,0,0,0 ,0,1,1,1,1,2
p2 := 2,3,3,3,3,3,3,3,3,3, 3,3,2,2,2,1,1,0,0,0 ,0,0,0,0,0,0,0,0,0,0 ,0,1,1,1,2,2
p3 := 2,3,3,3,3,3,3,3,3,3, 3,3,2,2,2,1,1,1,1,0 ,0,0,0,0,0,0,0,0,0,0 ,0,1,1,2,2,2
p4 := 2,2,3,3,3,3,3,3,3,3, 3,3,3,2,2,2,2,1,1,0 ,0,0,0,0,0,0,0,0,0,0 ,0,1,1,1,2,2
p5 := 2,2,3,3,3,3,3,3,3,3, 3,3,3,3,2,2,1,1,1,0 ,0,0,0,0,0,0,0,0,0,0 ,0,1,1,1,1,2
p6 := 2,2,3,3,3,3,3,3,3,3, 3,3,3,2,2,2,1,1,1,0 ,0,0,0,0,0,0,0,0,0,0 ,0,0,1,1,1,2
p7 := 2,2,2,3,3,3,3,3,3,3, 3,3,2,2,2,2,1,1,1,0 ,0,0,1,1,1,1,1,1,0,0 ,0,0,0,1,1,2
p8 := 2,2,2,2,3,3,3,3,3,3, 3,2,2,2,2,2,2,1,1,1 ,1,1,1,1,1,1,1,1,0,0 ,0,0,0,1,1,2
p9 := 2,2,2,2,2,3,2,2,3,2, 2,2,2,2,2,2,2,1,1,1 ,1,1,1,1,2,2,1,1,0,0 ,0,0,0,1,1,2
p10:= 2,2,2,2,2,2,2,2,2,2, 2,2,2,2,1,1,1,1,1,1 ,1,2,2,2,2,1,1,1,0,0 ,0,0,0,1,1,1
p11:= 2,2,2,2,2,2,2,2,2,2, 2,2,1,1,1,1,1,1,1,1 ,1,2,2,2,2,1,1,1,0,0 ,0,0,0,1,1,1
p12:= 1,2,2,2,2,2,2,2,2,2, 2,2,2,1,1,0,0,0,0,1 ,1,1,2,2,2,1,1,1,1,0 ,0,0,0,1,1,1
p13:= 1,1,2,2,2,2,2,2,2,2, 2,2,2,1,1,0,0,0,0,1 ,1,2,2,2,2,2,2,1,1,0 ,0,0,0,1,1,1
p14:= 1,1,2,2,2,2,2,2,2,2, 2,2,1,1,1,0,0,0,0,1 ,1,2,2,3,2,2,1,1,1,0 ,0,0,0,1,1,1
p15:= 1,2,2,2,2,2,2,2,2,2, 2,2,1,1,1,0,0,0,0,1 ,1,1,2,3,3,2,1,1,1,0 ,0,0,0,1,1,1
p16:= 1,2,2,2,2,2,2,3,2,2, 2,2,2,1,1,1,0,0,0,1 ,1,1,2,2,2,2,2,1,1,0 ,0,0,0,1,1,1
p17:= 1,2,2,2,2,2,3,3,3,3, 2,2,2,1,1,1,1,0,0,0 ,1,1,1,2,2,2,1,1,1,0 ,0,0,0,1,1,1
p18:= 1,2,2,2,2,2,3,3,3,3, 3,2,2,2,1,1,1,0,0,0 ,1,1,1,1,2,2,2,1,1,0 ,0,0,0,1,1,1
p19:= 1,2,2,2,2,3,3,3,3,3, 3,3,2,2,2,1,1,1,0,0 ,1,1,1,2,1,2,2,1,1,0 ,0,0,0,1,1,1
p20:= 1,2,2,2,2,2,3,3,3,3, 3,3,3,2,2,1,1,1,1,0 ,1,1,1,1,1,1,1,1,1,0 ,0,0,0,1,1,1
p21:= 1,1,2,2,2,2,3,3,3,3, 3,3,3,3,2,2,1,1,1,0 ,1,1,1,1,1,1,1,1,1,0 ,0,0,0,1,1,1
p22:= 1,1,2,2,2,2,2,3,3,3, 3,3,3,3,2,2,2,1,1,1 ,0,0,0,0,0,0,0,0,0,0 ,0,0,0,1,1,1
p23:= 1,1,1,2,2,2,2,2,3,3, 3,3,3,3,3,2,2,1,1,1 ,0,0,0,0,0,0,0,0,0,0 ,0,0,0,1,1,1
p24:= 1,2,2,2,2,2,2,2,2,3, 3,3,3,3,3,2,2,2,1,1 ,1,0,0,0,0,0,0,0,0,0 ,0,0,0,1,1,1
p25:= 2,2,1,2,2,2,2,2,2,2, 2,3,3,3,3,2,2,2,1,1 ,1,0,0,0,0,0,0,0,0,0 ,0,0,0,1,1,1
p26:= 1,1,1,2,2,2,2,2,2,2, 2,2,2,2,2,2,2,2,2,1 ,1,0,0,0,0,0,0,0,0,0 ,0,0,0,1,1,1
p27:= 1,1,1,1,1,2,2,2,2,2, 2,2,2,2,2,2,2,2,1,1 ,1,0,0,0,0,0,0,0,0,0 ,0,0,0,1,1,1
p28:= 0,1,1,1,1,2,2,2,2,2, 2,2,2,2,2,1,1,1,1,1 ,1,0,0,0,0,0,0,0,0,0 ,0,0,0,0,0,0
p29:= 0,0,0,1,1,1,1,1,1,1, 1,2,1,1,1,1,1,1,1,1 ,0,0,0,0,0,0,0,0,0,0 ,0,0,0,0,0,0
p30:= 0,0,0,1,1,1,1,2,2,1, 1,1,1,1,1,1,1,0,0,0 ,0,0,0,0,0,0,0,0,0,0 ,0,0,0,0,0,0
p31:= 1,1,1,1,1,2,2,1,1,1, 1,1,1,1,0,0,0,0,0,0 ,0,1,1,1,1,1,0,0,0,0 ,0,0,0,0,0,0
p32:= 1,1,1,1,1,1,2,2,2,2, 1,1,1,1,0,0,0,0,0,0 ,0,1,1,1,1,1,0,0,0,0 ,0,0,0,0,0,0
p33:= 1,2,2,2,1,2,2,2,2,2, 1,1,1,1,1,0,0,0,0,0 ,0,1,1,2,1,1,1,0,0,0 ,0,0,0,0,0,0
p34:= 1,2,2,2,2,2,2,2,2,2, 2,1,1,1,1,0,0,0,0,0 ,0,1,1,2,1,1,1,0,0,0 ,0,0,0,1,1,1
p35:= 1,2,2,2,2,2,2,2,2,2, 2,2,2,1,1,0,0,0,0,0 ,0,1,1,2,2,1,1,0,0,0 ,0,0,0,1,1,1
p36:= 2,2,2,2,2,2,2,2,2,2, 2,2,2,1,1,0,0,0,0,0 ,0,1,1,1,1,1,1,0,0,0 ,0,0,0,1,1,2

aas:='LEU GLY ALA VAL ILE PRO SER THR GLU ASP GLN ASN LYS ARG HIS PHE TYR TRP MET CYS'

print "X0=0 Y0=-36 X1=36 Y1=0" >$file(1:indexr('$file','.'))grf
print "x0=-200 y0=-50 x1=200 y1=350" >>
print "border=off font=Times textsize=18"  >>

# ------------------------------------------------------------- draw background

if (.not.nobackground) then
  print "color=blue"  >>
  fillp:=5,4,2
  do y 1 36
    lastc=p{$y}(1)
    lastx=1
    do x 1 36
      i=p{$y}($x)
      if (i.ne.lastc) then
        if (lastc.gt.0) then
          print "fill=$fillp(lastc)" >>
          print "rectangle ${lastx-1} ${1-y} ${x-1} ${-y}" >>
        end if
        lastx=x
        lastc=i
      end if
    end do
    if (i.gt.0) then
      print "fill=$fillp(i)" >>
      print "rectangle ${lastx-1} ${1-y} ${x-1} ${-y}" >>
    end if
  end do
end if

# -------------------------------------------------------- draw angle positions

print "color=black textsize=18 autoscale=off" >>
print "X0=-180 Y0=-180 X1=180 Y1=180" >>
print "label x -180 -180" >>
print "label x -90 -90"   >>
print "label x 0 0"       >>
print "label x 90 90"     >>
print "label x 180 180"   >>
print "label y -180 -180" >>
print "label y -90 -90"   >>
print "label y 0 0"       >>
print "label y 90 90"     >>
print "label y 180 180"   >>
print "text -220 -5 @Sy"  >>
print "frame @Sf tics=off" >>

print "    Struct    fav    add    gen    dis"
print "    ------    ---    ---    ---    ---"

s0=0; s1=0; s2=0; s3=0; n=0
do k 1 nstruct
  struct=istruct(k)
  structure copy struct 0
  c0=0; c1=0; c2=0; c3=0; str:=
  do i 2 nr-1
    rp=rnam(i-1); r=rnam(i); rn=rnam(i+1)
    if (mod(index(aas,rp(1:3)),4).eq.1 .and. mod(index(aas,r(1:3)),4).eq.1 .and. \
        mod(index(aas,rn(1:3)),4).eq.1 .and. r.ne.'GLY' .and. r(1:3).ne.'PRO') then
      j=rnum(i)
      phi=dval('PHI $j')
      psi=dval('PSI $j')
      phi=mod(phi,360.0)
      if (abs(phi).ge.180.0) phi=phi-sign(360.0,phi)
      psi=mod(psi,360.0)
      if (abs(psi).ge.180.0) psi=psi-sign(360.0,psi)
      print "$phi $psi" >>
      xx=max(1,min(36,1+int((phi+180.0)/10)))
      yy=max(1,min(36,36-int((psi+180.0)/10)))
      x=p{$yy}($xx)
      if (x.eq.0) then
        c0=c0+1
        if (label) print "text $phi ${psi+6} $r $j" >>
        if (c0.eq.1) then
          str:="($r $j"
        else if (c0.le.10) then
          str:="$str, $r $j"
        else if (c0.eq.11) then
          str:="$str,..."
        end if
      else if (x.eq.1) then
        c1=c1+1
      else if (x.eq.2) then
        c2=c2+1
      else
        c3=c3+1
      end if
      n=n+1
    end if
  end do
  print "border=on textsize=10 marksize=4 fill=4" >>
  print "dot x y1" >>
  if (c0.gt.0) str:="$str)"
  print "    $struct(I4)  $c3(I7)$c2(I7)$c1(I7)$c0(I7)  $str"
  s0=s0+c0; s1=s1+c1; s2=s2+c2; s3=s3+c3
end do
n=max(1,n);
s0=100.0*s0/n; s1=100.0*s1/n; s2=100.0*s2/n; s3=100.0*s3/n
print "     all    $s3(F6.1)\%$s2(F6.1)\%$s1(F6.1)\%$s0(F6.1)\%"
print "textsize=18"  >>
print "clip off"  >>
print "text -200 -280 $s3(F6.1)\% in most favored regions" >>
print "text -200 -300 $s2(F6.1)\% in additionally allowed regions" >>
print "text -200 -320 $s1(F6.1)\% in generously allowed regions" >>
print "text -200 -340 $s0(F6.1)\% in disallowed regions" >>
print " " >>.

graf $file replace