! Module file: acceptavesubave
!
! CNS MODULE
! **********
!
! Authors: Gregory L. Warren and Axel T. Brunger
!
! copyright Yale University
!
! version 03/05/98
!
! Function:
!          Calculate average values for accepted or trial
!          ensembles
!
! Requirements:
!       Use symbols initilized and used in the macros initave 
!       and printaccept.             
!
!

module {acceptavesubave}
(
   &ave=ave;                 {INPUT/OUPUT: summed value from acceptsubave}
   &ave2=ave2;               {INPUT/OUPUT: squared sum values from acceptsubave}
   &count=count;             {INPUT: number of accepted or trial structures}
   &nmr=nmr;                 {INPUT: parameters for nmr restraint data}
)

checkversion 1.3

set message ? end
evaluate ($message_old=$result)
set echo ? end
evaluate ($echo_old=$result)
if ( $log_level = verbose ) then
  set echo=on message=normal end
else
  set echo=off message=off end
end if

   if (&count > 1) then
      {- rms sd values -}
      evaluate (&ave2.rms.bond=sqrt
               ((&ave2.rms.bond-(&ave.rms.bond^2/&count))
               /(&count-1)))
      evaluate (&ave2.rms.angl=sqrt
               ((&ave2.rms.angl-(&ave.rms.angl^2/&count))
               /(&count-1)))
      evaluate (&ave2.rms.impr=sqrt
               ((&ave2.rms.impr-(&ave.rms.impr^2/&count))
               /(&count-1)))
      evaluate (&ave2.rms.dihe=sqrt
               ((&ave2.rms.dihe-(&ave.rms.dihe^2/&count))
               /(&count-1)))
      evaluate (&ave2.rms.noe=sqrt
               ((&ave2.rms.noe-(&ave.rms.noe^2/&count))
               /(&count-1)))
      evaluate (&ave2.rms.jcoup=sqrt 
               ((&ave2.rms.jcoup-(&ave.rms.jcoup^2/&count)) 
               /(&count-1)))
      evaluate (&ave2.rms.oneb=sqrt 
               ((&ave2.rms.oneb-(&ave.rms.oneb^2/&count)) 
               /(&count-1)))
      evaluate (&ave2.rms.carb.a=sqrt
               ((&ave2.rms.carb.a-(&ave.rms.carb.a^2/&count))
               /(&count-1)))
      evaluate (&ave2.rms.carb.b=sqrt
               ((&ave2.rms.carb.b-(&ave.rms.carb.b^2/&count))
               /(&count-1)))
      evaluate (&ave2.rms.prot.all=sqrt
               ((&ave2.rms.prot.all-(&ave.rms.prot.all^2/&count))
               /(&count-1)))
      evaluate (&ave2.rms.dani=sqrt
               ((&ave2.rms.dani-(&ave.rms.dani^2/&count))
               /(&count-1)))
      evaluate (&ave2.rms.sani=sqrt
               ((&ave2.rms.sani-(&ave.rms.sani^2/&count))
               /(&count-1)))
      evaluate (&ave2.rms.cdih=sqrt
               ((&ave2.rms.cdih-(&ave.rms.cdih^2/&count))
               /(&count-1)))

      {- violation sd values -}
      evaluate (&ave2.viol.bond=sqrt
               ((&ave2.viol.bond-(&ave.viol.bond^2/&count))
               /(&count-1)))
      evaluate (&ave2.viol.angl=sqrt
               ((&ave2.viol.angl-(&ave.viol.angl^2/&count))
               /(&count-1)))
      evaluate (&ave2.viol.impr=sqrt
               ((&ave2.viol.impr-(&ave.viol.impr^2/&count))
               /(&count-1)))
      evaluate (&ave2.viol.dihe=sqrt
               ((&ave2.viol.dihe-(&ave.viol.dihe^2/&count))
               /(&count-1)))
      evaluate (&ave2.viol.vdw=sqrt
               ((&ave2.viol.vdw-(&ave.viol.vdw^2/&count))
               /(&count-1)))
      evaluate (&ave2.viol.noe=sqrt
               ((&ave2.viol.noe-(&ave.viol.noe^2/&count))
               /(&count-1)))
      evaluate (&ave2.viol.jcoup=sqrt 
               ((&ave2.viol.jcoup-(&ave.viol.jcoup^2/&count)) 
               /(&count-1)))
      evaluate (&ave2.viol.oneb=sqrt 
               ((&ave2.viol.oneb-(&ave.viol.oneb^2/&count)) 
               /(&count-1)))
      evaluate (&ave2.viol.carb=sqrt
               ((&ave2.viol.carb-(&ave.viol.carb^2/&count))
               /(&count-1)))
      evaluate (&ave2.viol.prot.all=sqrt
               ((&ave2.viol.prot.all-(&ave.viol.prot.all^2/&count))
               /(&count-1)))
      evaluate (&ave2.viol.dani=sqrt
               ((&ave2.viol.dani-(&ave.viol.dani^2/&count))
               /(&count-1)))
      evaluate (&ave2.viol.sani=sqrt
               ((&ave2.viol.sani-(&ave.viol.sani^2/&count))
               /(&count-1)))
      evaluate (&ave2.viol.cdih=sqrt
               ((&ave2.viol.cdih-(&ave.viol.cdih^2/&count))
               /(&count-1)))

      evaluate (&ave2.rotx.plan=sqrt
               ((&ave2.rotx.plan-(&ave.rotx.plan^2/&count))
               /(&count-1)))
      evaluate (&ave2.roty.plan=sqrt
               ((&ave2.roty.plan-(&ave.roty.plan^2/&count))
               /(&count-1)))
      evaluate (&ave2.rotz.plan=sqrt
               ((&ave2.rotz.plan-(&ave.rotz.plan^2/&count))
               /(&count-1)))

      {- proton chemical shift class rms and viol sd values -}
      evaluate ($count=1)
      while (&exist%nmr.prot.file.$count=true) loop nloop
         evaluate ($clsname = "P"+encode($count))
         if (&nmr.prot.file.$count # "") then
            evaluate (&ave2.rms.prot.$clsname=sqrt
                     ((&ave2.rms.prot.$clsname-(&ave.rms.prot.$clsname^2/&count))
                     /(&count-1)))
            evaluate (&ave2.viol.prot.$clsname=sqrt
                     ((&ave2.viol.prot.$clsname-(&ave.viol.prot.$clsname^2/&count))
                     /(&count-1)))
         end if
         evaluate ($count = $count + 1)
      end loop nloop 
   else
      {- rms sd values -}
      evaluate (&ave2.rms.bond=0.)
      evaluate (&ave2.rms.angl=0.)
      evaluate (&ave2.rms.impr=0.)
      evaluate (&ave2.rms.dihe=0.)
      evaluate (&ave2.rms.noe=0.)
      evaluate (&ave2.rms.jcoup=0.)
      evaluate (&ave2.rms.oneb=0.)
      evaluate (&ave2.rms.carb.a=0.)
      evaluate (&ave2.rms.carb.b=0.)
      evaluate (&ave2.rms.prot.all=0.)
      evaluate (&ave2.rms.dani=0.)
      evaluate (&ave2.rms.sani=0.)
      evaluate (&ave2.rms.cdih=0.)

      {- violation sd values -}
      evaluate (&ave2.viol.bond=0.)
      evaluate (&ave2.viol.angl=0.)
      evaluate (&ave2.viol.impr=0.)
      evaluate (&ave2.viol.dihe=0.)
      evaluate (&ave2.viol.vdw=0.)
      evaluate (&ave2.viol.noe=0.)
      evaluate (&ave2.viol.jcoup=0.)
      evaluate (&ave2.viol.oneb=0.)
      evaluate (&ave2.viol.carb=0.)
      evaluate (&ave2.viol.prot.all=0.)
      evaluate (&ave2.viol.dani=0.)
      evaluate (&ave2.viol.sani=0.)
      evaluate (&ave2.viol.cdih=0.)

      evaluate (&ave2.rotx.plan=0.)
      evaluate (&ave2.roty.plan=0.)
      evaluate (&ave2.rotz.plan=0.)

      {- proton chemical shift class rms and viol sd values -}
      evaluate ($count=1)
      while (&exist%nmr.prot.file.$count=true) loop nloop
         evaluate ($clsname = "P"+encode($count))
         if (&nmr.prot.file.$count # "") then
            evaluate (&ave2.rms.prot.$clsname=0.)
            evaluate (&ave2.viol.prot.$clsname=0.)
         end if
         evaluate ($count = $count + 1)
      end loop nloop 
   end if      

   {- average rms values -}
      evaluate (&ave.rms.bond=&ave.rms.bond/&count)
      evaluate (&ave.rms.angl=&ave.rms.angl/&count)
      evaluate (&ave.rms.impr=&ave.rms.impr/&count)
      evaluate (&ave.rms.dihe=&ave.rms.dihe/&count)
      evaluate (&ave.rms.noe=&ave.rms.noe/&count)
      evaluate (&ave.rms.jcoup=&ave.rms.jcoup/&count) 
      evaluate (&ave.rms.oneb=&ave.rms.oneb/&count) 
      evaluate (&ave.rms.carb.a=&ave.rms.carb.a/&count)
      evaluate (&ave.rms.carb.b=&ave.rms.carb.b/&count)
      evaluate (&ave.rms.prot.all=&ave.rms.prot.all/&count)
      evaluate (&ave.rms.dani=&ave.rms.dani/&count) 
      evaluate (&ave.rms.sani=&ave.rms.sani/&count) 
      evaluate (&ave.rms.cdih=&ave.rms.cdih/&count) 

   {- average violation values -}
      evaluate (&ave.viol.bond=&ave.viol.bond/&count)
      evaluate (&ave.viol.angl=&ave.viol.angl/&count)
      evaluate (&ave.viol.impr=&ave.viol.impr/&count)
      evaluate (&ave.viol.dihe=&ave.viol.dihe/&count)
      evaluate (&ave.viol.vdw=&ave.viol.vdw/&count)
      evaluate (&ave.viol.noe=&ave.viol.noe/&count)
      evaluate (&ave.viol.jcoup=&ave.viol.jcoup/&count) 
      evaluate (&ave.viol.oneb=&ave.viol.oneb/&count) 
      evaluate (&ave.viol.carb=&ave.viol.carb/&count)
      evaluate (&ave.viol.prot.all=&ave.viol.prot.all/&count)
      evaluate (&ave.viol.dani=&ave.viol.dani/&count) 
      evaluate (&ave.viol.sani=&ave.viol.sani/&count) 
      evaluate (&ave.viol.cdih=&ave.viol.cdih/&count) 

      evaluate (&ave.rotx.plan=&ave.rotx.plan/&count) 
      evaluate (&ave.roty.plan=&ave.roty.plan/&count) 
      evaluate (&ave.rotz.plan=&ave.rotz.plan/&count) 

   {- proton chemical shift class rms and viol average values -}
      evaluate ($count=1)
      while (&exist%nmr.prot.file.$count=true) loop nloop
         evaluate ($clsname = "P"+encode($count))
         if (&nmr.prot.file.$count # "") then
            evaluate (&ave.rms.prot.$clsname=&ave.rms.prot.$clsname/&count)
            evaluate (&ave.viol.prot.$clsname=&ave.viol.prot.$clsname/&count)
         end if
         evaluate ($count = $count + 1)
      end loop nloop

set message=$message_old echo=$echo_old end