remarks Analysis protocol

set message on echo on end

{*==========================================================================*}
{*====================== SET FILENAMES AND VARIABLES =======================*}
{*==========================================================================*}

{* type of non-bonded parameters *}
{* choice: "PROLSQ" "PARMALLH6" "PARALLHDG" "OPLSX" *}
{* The water refinement uses the OPLSX parameters *}

evaluate ( $par_nonbonded = "OPLSX" )

{* parameter file(s) *}
evaluate ( $par.1 = "TOPPAR:parallhdg5.3.pro" )
evaluate ( $par.2 = "TOPPAR:parallhdg5.3.sol" )
evaluate ( $par.3 = "" )
evaluate ( $par.4 = "" )
evaluate ( $par.5 = "" )

{* solvent topology file *}
evaluate ( $solvent_topology = "TOPPAR:topallhdg5.3.sol" )

{* structure file(s) *}
evaluate ( $struct.1 = "STRUCTURES:1xxx.psf" )
evaluate ( $struct.2 = "" )
evaluate ( $struct.3 = "" )
evaluate ( $struct.4 = "" )
evaluate ( $struct.5 = "" )

{* input coordinate file(s) *}
evaluate ( $pdb.in.file.1 = "INPUTCOORDINATES:1xxx_13.pdb" )

{* output coordinate file(s) *}
evaluate ( $pdb.out.file.1 = "OUTPUTCOORDINATES:analyzed_13.pdb" )

{* input distance restraint file(s) *}
evaluate ( $noe.file.1 = "TABLES:1xxx_noe.tbl" )

{* Averaging for NOE restraints *}
evaluate ( $noe.ave = sum )

{* input dihedral restraint file(s) *}
evaluate ( $cdih.file.1 = "TABLES:1xxx_dihe.tbl" )

{* seed for velocity generation *}
evaluate ( $seed = 1234123 )

{* set number of md steps for the heating stage *}
evaluate ( $mdsteps.heat = 20 )

{* set number of md steps for the hot md stage *}
evaluate ( $mdsteps.hot = 100 )

{* set number of md steps for the cooling stage *}
evaluate ( $mdsteps.cool = 20 )

{* seed for velocity generation *}
evaluate ( $seed = 1234123 )

{* set acceptance criteria *}
evaluate ( $accept.noe = 0.50 )
evaluate ( $accept.cdih = 5.00 )
evaluate ( $accept.coup = 1.00 )
evaluate ( $accept.sani = 0.00 )
evaluate ( $accept.vean = 5.00 )


{*==========================================================================*}
{*=========== READ THE STRUCTURE, TOPOLOGY AND PARAMETER FILES =============*}
{*==========================================================================*}

structure @$struct.1 end
topology  @$solvent_topology end
parameter @$par.1 @$par.2 end

{*==========================================================================*}
{*================== SET VALUES FOR NONBONDED PARAMETERS ===================*}
{*==========================================================================*}

parameter
  nbonds
    nbxmod=5 atom cdiel shift
    cutnb=9.5 ctofnb=8.5 ctonnb=6.5 eps=1.0 e14fac=0.4 inhibit 0.25
    wmin=0.5
    tolerance  0.5
  end
end


{*==========================================================================*}
{*========================= READ THE COORDINATES ===========================*}
{*==========================================================================*}

coor @@$pdb.in.file.1

{*==========================================================================*}
{*========================= READ THE EXPERIMENTAL DATA =====================*}
{*==========================================================================*}

noe reset
  nrestraints = 10000
  ceiling     = 100

  @$noe.file.1
  averaging  * $noe.ave
  potential  * soft
  scale      * 50
  sqconstant * 1.0
  sqexponent * 2
  soexponent * 1
  rswitch    * 1.0
  sqoffset   * 0.0
  asymptote  * 2.0
end

restraints dihedral reset
   @$cdih.file.1
   scale = 200
end

{*==========================================================================*}
{*============================ SET FLAGS ===================================*}
{*==========================================================================*}

flags exclude *
      include bond angle dihe impr vdw elec
              noe cdih coup oneb carb ncs dani
              vean sani prot harm
end

  constraints
    interaction (not resname ANI) (not resname ANI) weights * 1 dihed 5 end
    interaction ( resname ANI) ( resname ANI) weights * 1  end
  end

{*==========================================================================*}
{*======================= CHECK RESTRAINT VIOLATIONS =======================*}
{*==========================================================================*}

evaluate ($accept = 0)

print threshold=$accept.noe noe
evaluate ($v_noe = $violations)
print threshold=0.5 noe
evaluate ($v_noe_0.5 = $violations)
print threshold=0.4 noe
evaluate ($v_noe_0.4 = $violations)
print threshold=0.3 noe
evaluate ($v_noe_0.3 = $violations)
print threshold=0.2 noe
evaluate ($v_noe_0.2 = $violations)
print threshold=0.1 noe
evaluate ($v_noe_0.1 = $violations)
evaluate ($rms_noe = $result)

print threshold=$accept.cdih cdih
evaluate ($rms_cdih=$result)
evaluate ($v_cdih = $violations)

coupl print thres = $accept.coup class c1 end
evaluate ($rms_coup = $result)
evaluate ($v_coup = $violations)

sani print threshold = $accept.sani class rdc1 end
evaluate ($rms_sani = $result)
evaluate ($v_sani = $violations)

vean print threshold = $accept.vean class vea1 end
evaluate( $rms_vean = $result)
evaluate( $v_vean = $violations)

if ($v_noe > 0) then  evaluate ( $accept=$accept + 1 ) end if
if ($v_cdih > 0) then  evaluate ( $accept=$accept + 1 ) end if
if ($v_coup > 0) then  evaluate ( $accept=$accept + 1 ) end if
if ($v_sani > 0) then  evaluate ( $accept=$accept + 1 ) end if
if ($v_vean > 0) then  evaluate ( $accept=$accept + 1 ) end if

if ($accept = 0 ) then
  evaluate ( $label = "ACCEPTED" )
else
  evaluate ( $label = "NOT ACCEPTED" )
end if

{*==========================================================================*}
{*================ COMPUTE RMS DIFFERENCE BETWEEN OBSERVED =================*}
{*========================= AND MODEL DISTANCES. ===========================*}
{*==========================================================================*}

constraints interaction
   (not resname TIP* and not resname ANI)
   (not resname TIP* and not resname ANI)
end

print thres=0.05 bonds
evaluate ($rms_bonds=$result)
print thres=5. angles
evaluate ($rms_angles=$result)
print thres=5. impropers
evaluate ($rms_impropers=$result)

energy end

remarks Structure $label
remarks E-overall:                             $ener
remarks E-NOE_restraints:                      $noe
remarks E-CDIH_restraints:                     $cdih
remarks E-COUP_restraints:                     $coup
remarks E-SANI_restraints:                     $sani
remarks E-VEAN_restraints:                     $vean
remarks RMS-NOE_restraints:                    $rms_noe
remarks RMS-CDIH_restraints:                   $rms_cdih
remarks RMS-COUP_restraints:                   $rms_coup
remarks RMS-SANI_restraints:                   $rms_sani
remarks RMS-VEAN_restraints:                   $rms_vean
remarks NOE Acceptance criterium:              $accept.noe
remarks NOE violations > $accept.noe:          $v_noe
remarks All NOE viol. (>0.5,0.4,0.3,0.2,0.1A): $v_noe_0.5 $v_noe_0.4 $v_noe_0.3 $v_noe_0.2 $v_noe_0.1
remarks CDIH Acceptance criterium:             $accept.cdih
remarks CDIH violations > $accept.cdih:        $v_cdih
remarks COUP Acceptance criterium:             $accept.coup
remarks COUP violations > $accept.coup:        $v_coup
remarks SANI Acceptance criterium:             $accept.sani
remarks SANI violations > $accept.sani:        $v_sani
remarks VEAN Acceptance criterium:             $accept.vean
remarks VEAN violations > $accept.vean:        $v_vean

write coordinates sele= (not resn TIP3) output =$pdb.out.file.1 end

stop