{+ file: dna-rna_restraints.def  +}
{+ directory: auxiliary +}
{+ description: Creates restraints to maintain conformation of DNA/RNA +}
{+ comment: 
           This file is to be read by refinement files that modify
           atom coordinates +}
{+ authors: Axel T. Brunger, and Paul D. Adams +}
{+ copyright: Yale University +}

{- begin block parameter definition -} define(

{========================== base planarity ===========================}

{* Nucleic acid residues to have base planarity restrained. This selection
   must only include nucleotide residues *}
{===>} bases_planar=(resid 402:419 or resid 421:438);

{========================== sugar puckers ============================}

{* residues with sugar pucker restrained - group 1 *}
{===>} pucker_1=(resid 402:419);

{* sugar pucker conformation of group 1 *}
{+ choice: "a-form" "b-form" "other" +}
{===>} form_1="a-form";

{* user defined sugar pucker for group 1 *}

{* dihedral C1'-C2'-C3'-C4' *}
{===>} dihedral_nu2_1=37.053;
{* dihedral C5'-C4'-C3'-C2' *}
{===>} dihedral_nu3_1=-155.59;
{* dihedral C1'-O4'-C4'-C5' *}
{===>} dihedral_nu4_1=144.26;

{* residues with sugar pucker restrained - group 2 *}
{===>} pucker_2=(resid 421:438);

{* sugar pucker conformation of group 2 *}
{+ choice: "a-form" "b-form" "other" +}
{===>} form_2="b-form";

{* user defined sugar pucker for group 2 *}

{* dihedral C1'-C2'-C3'-C4' *}
{===>} dihedral_nu2_2=37.053;
{* dihedral C5'-C4'-C3'-C2' *}
{===>} dihedral_nu3_2=-155.59;
{* dihedral C1'-O4'-C4'-C5' *}
{===>} dihedral_nu4_2=144.26;

{* residues with sugar pucker restrained - group 3 *}
{===>} pucker_3=(none);

{* sugar pucker conformation of group 3 *}
{+ choice: "a-form" "b-form" "other" +}
{===>} form_3="a-form";

{* user defined sugar pucker for group 3 *}

{* dihedral C1'-C2'-C3'-C4' *}
{===>} dihedral_nu2_3=37.053;
{* dihedral C5'-C4'-C3'-C2' *}
{===>} dihedral_nu3_3=-155.59;
{* dihedral C1'-O4'-C4'-C5' *}
{===>} dihedral_nu4_3=144.26;

{* residues with sugar pucker restrained - group 4 *}
{===>} pucker_4=(none);

{* sugar pucker conformation of group 4 *}
{+ choice: "a-form" "b-form" "other" +}
{===>} form_4="a-form";

{* user defined sugar pucker for group 4 *}

{* dihedral C1'-C2'-C3'-C4' *}
{===>} dihedral_nu2_4=37.053;
{* dihedral C5'-C4'-C3'-C2' *}
{===>} dihedral_nu3_4=-155.59;
{* dihedral C1'-O4'-C4'-C5' *}
{===>} dihedral_nu4_4=144.26;

{====================== Watson-Crick base pairs ======================}

{* residues which form Watson-Crick pairs *}

{* selection for pair 1 base A *}
{===>} base_a_1=(resid 402);
{* selection for pair 1 base B *}
{===>} base_b_1=(resid 438);

{* selection for pair 2 base A *}
{===>} base_a_2=(resid 403);
{* selection for pair 2 base B *}
{===>} base_b_2=(resid 437);

{* selection for pair 3 base A *}
{===>} base_a_3=(resid 404);
{* selection for pair 3 base B *}
{===>} base_b_3=(resid 436);

{* selection for pair 4 base A *}
{===>} base_a_4=(resid 405);
{* selection for pair 4 base B *}
{===>} base_b_4=(resid 435);

{* selection for pair 5 base A *}
{===>} base_a_5=(resid 406);
{* selection for pair 5 base B *}
{===>} base_b_5=(resid 434);

{* selection for pair 6 base A *}
{===>} base_a_6=(resid 407);
{* selection for pair 6 base B *}
{===>} base_b_6=(resid 433);

{* selection for pair 7 base A *}
{===>} base_a_7=(resid 408);
{* selection for pair 7 base B *}
{===>} base_b_7=(resid 432);

{* selection for pair 8 base A *}
{===>} base_a_8=(resid 409);
{* selection for pair 8 base B *}
{===>} base_b_8=(resid 431);

{* selection for pair 9 base A *}
{===>} base_a_9=(resid 410);
{* selection for pair 9 base B *}
{===>} base_b_9=(resid 430);

{* selection for pair 10 base A *}
{===>} base_a_10=(resid 411);
{* selection for pair 10 base B *}
{===>} base_b_10=(resid 429);

{* selection for pair 11 base A *}
{===>} base_a_11=(resid 412);
{* selection for pair 11 base B *}
{===>} base_b_11=(resid 428);

{* selection for pair 12 base A *}
{===>} base_a_12=(resid 413);
{* selection for pair 12 base B *}
{===>} base_b_12=(resid 427);

{* selection for pair 13 base A *}
{===>} base_a_13=(resid 414);
{* selection for pair 13 base B *}
{===>} base_b_13=(resid 426);

{* selection for pair 14 base A *}
{===>} base_a_14=(resid 415);
{* selection for pair 14 base B *}
{===>} base_b_14=(resid 425);

{* selection for pair 15 base A *}
{===>} base_a_15=(resid 416);
{* selection for pair 15 base B *}
{===>} base_b_15=(resid 424);

{* selection for pair 16 base A *}
{===>} base_a_16=(resid 417);
{* selection for pair 16 base B *}
{===>} base_b_16=(resid 423);

{* selection for pair 17 base A *}
{===>} base_a_17=(resid 418);
{* selection for pair 17 base B *}
{===>} base_b_17=(resid 422);

{* selection for pair 18 base A *}
{===>} base_a_18=(resid 419);
{* selection for pair 18 base B *}
{===>} base_b_18=(resid 421);

{* selection for pair 19 base A *}
{===>} base_a_19=(none);
{* selection for pair 19 base B *}
{===>} base_b_19=(none);

{* selection for pair 20 base A *}
{===>} base_a_20=(none);
{* selection for pair 20 base B *}
{===>} base_b_20=(none);

{* selection for pair 21 base A *}
{===>} base_a_21=(none);
{* selection for pair 21 base B *}
{===>} base_b_21=(none);

{* selection for pair 22 base A *}
{===>} base_a_22=(none);
{* selection for pair 22 base B *}
{===>} base_b_22=(none);

{* selection for pair 23 base A *}
{===>} base_a_23=(none);
{* selection for pair 23 base B *}
{===>} base_b_23=(none);

{* selection for pair 24 base A *}
{===>} base_a_24=(none);
{* selection for pair 24 base B *}
{===>} base_b_24=(none);

{* selection for pair 25 base A *}
{===>} base_a_25=(none);
{* selection for pair 25 base B *}
{===>} base_b_25=(none);

{* selection for pair 26 base A *}
{===>} base_a_26=(none);
{* selection for pair 26 base B *}
{===>} base_b_26=(none);

{* selection for pair 27 base A *}
{===>} base_a_27=(none);
{* selection for pair 27 base B *}
{===>} base_b_27=(none);

{* selection for pair 28 base A *}
{===>} base_a_28=(none);
{* selection for pair 28 base B *}
{===>} base_b_28=(none);

{* selection for pair 29 base A *}
{===>} base_a_29=(none);
{* selection for pair 29 base B *}
{===>} base_b_29=(none);

{* selection for pair 30 base A *}
{===>} base_a_30=(none);
{* selection for pair 30 base B *}
{===>} base_b_30=(none);

{* selection for pair 31 base A *}
{===>} base_a_31=(none);
{* selection for pair 31 base B *}
{===>} base_b_31=(none);

{* selection for pair 32 base A *}
{===>} base_a_32=(none);
{* selection for pair 32 base B *}
{===>} base_b_32=(none);

{===========================================================================}
{        things below this line do not normally need to be changed          }
{===========================================================================}

 ) {- end block parameter definition -}

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

 checkversion 1.2

 {- the planarity restraints -}

 for $id in id ( &bases_planar and tag ) loop plan

   show (segid) (id $id)
   evaluate ($segid=$result)
   show (resid) (id $id)
   evaluate ($resid=$result)

   restraints plane

     group
       selection=( segid $segid and resid $resid and 
                  (resname THY or resname CYT or resname GUA or
                   resname ADE or resname URI) and
                   not ( name c#' or name h#' or name o#p or
                         name o#' or name p or name h#t or name o#t ) )
       weight=400.0
     end

   end

 end loop plan

 flags include plan end

 {- the dihedral restraints for sugar puckers -}

 restraints dihedral          
   nassign=1000
 end

 evaluate ($group=1)
 evaluate ($done=false)

 while ( $done = false ) loop dihe

   if ( &exist_pucker_$group = true ) then

     show sum(1) ( &pucker_$group )

     if ( $result > 0 ) then

       if ( &form_$group = "a-form" ) then
         evaluate ($dihedral_1=37.053)
         evaluate ($dihedral_2=-155.59)
         evaluate ($dihedral_3=144.26)
       elseif ( &form_$group = "b-form" ) then
         evaluate ($dihedral_1=-34.9)
         evaluate ($dihedral_2=-86.4)
         evaluate ($dihedral_3=106.4)
       elseif ( &form_$group = "other" ) then
         evaluate ($dihedral_1=&dihedral_nu2_$group)
         evaluate ($dihedral_2=&dihedral_nu3_$group)
         evaluate ($dihedral_3=&dihedral_nu4_$group)
       end if

       for $id in id ( &pucker_$group and tag ) loop resid

         show (segid) (id $id)
         evaluate ($segid=$result)
         show (resid) ( id $id )
         evaluate ($resid=$result)

         restraints dihedral

           assign  ( segid $segid and resid $resid and name c1' )
                   ( segid $segid and resid $resid and name c2' )
                   ( segid $segid and resid $resid and name c3' )
                   ( segid $segid and resid $resid and name c4' ) 
                                                       20.0 $dihedral_1 0.0 2

           assign  ( segid $segid and resid $resid and name c5' )
                   ( segid $segid and resid $resid and name c4' )
                   ( segid $segid and resid $resid and name c3' )
                   ( segid $segid and resid $resid and name c2' ) 
                                                       20.0 $dihedral_2 0.0 2

           assign  ( segid $segid and resid $resid and name c1' )
                   ( segid $segid and resid $resid and name o4' )
                   ( segid $segid and resid $resid and name c4' )
                   ( segid $segid and resid $resid and name c5' ) 
                                                       20.0 $dihedral_3 0.0 2
  
           scale=20.0

         end

       end loop resid

     end if

   else

     evaluate ($done=true)

   end if

   evaluate ($group=$group+1)

 end loop dihe

 flags include cdih end

 {- Watson-Crick base pairing -}

 noe
   nres=5000
   class hres
   ceiling=1000
   averaging hres cent
   potential hres square
   sqconstant hres 1.
   sqexponent hres 2
   scale hres 50.
 end           

 evaluate ($pair=1)
 evaluate ($done=false)

 while ( $done = false ) loop noe

   if ( &exist_base_a_$pair = true ) then

     if ( &exist_base_b_$pair = true ) then

       show ( resname ) ( &base_a_$pair and name C1' ) 
       evaluate ($ares=$result)
       show ( resname ) ( &base_b_$pair and name C1' ) 
       evaluate ($bres=$result)

       noe

         if ( $ares = THY ) then

           assign (&base_a_$pair and name o4) 
                  (&base_b_$pair and name n6) 2.95 0.01 0.01 
           assign (&base_a_$pair and name n3) 
                  (&base_b_$pair and name n1) 2.82 0.01 0.01
           assign (&base_a_$pair and name o4) 
                  (&base_b_$pair and name n1) 3.63 0.01 0.01 
           assign (&base_a_$pair and name o2) 
                  (&base_b_$pair and name n6) 5.40 0.01 0.01

         elseif ( $ares = URI ) then

           assign (&base_a_$pair and name o4) 
                  (&base_b_$pair and name n6) 2.95 0.01 0.01 
           assign (&base_a_$pair and name n3) 
                  (&base_b_$pair and name n1) 2.82 0.01 0.01
           assign (&base_a_$pair and name o4) 
                  (&base_b_$pair and name n1) 3.63 0.01 0.01 
           assign (&base_a_$pair and name o2) 
                  (&base_b_$pair and name n6) 5.40 0.01 0.01

         elseif ( $ares = ADE ) then

             assign (&base_b_$pair and name o4) 
                    (&base_a_$pair and name n6) 2.95 0.01 0.01 
             assign (&base_b_$pair and name n3) 
                    (&base_a_$pair and name n1) 2.82 0.01 0.01
             assign (&base_b_$pair and name o4) 
                    (&base_a_$pair and name n1) 3.63 0.01 0.01 
             assign (&base_b_$pair and name o2) 
                    (&base_a_$pair and name n6) 5.40 0.01 0.01

         elseif ( $ares = CYT ) then

           assign (&base_a_$pair and name n4) 
                  (&base_b_$pair and name o6) 2.91 0.01 0.01 
           assign (&base_a_$pair and name n3) 
                  (&base_b_$pair and name n1) 2.95 0.01 0.01
           assign (&base_a_$pair and name o2)
                  (&base_b_$pair and name n2) 2.86 0.01 0.01
           assign (&base_a_$pair and name n3)
                  (&base_b_$pair and name n2) 3.65 0.01 0.01 
           assign (&base_a_$pair and name o2)
                  (&base_b_$pair and name o6) 5.42 0.01 0.01

         elseif ( $ares = GUA ) then

           assign (&base_b_$pair and name n4) 
                  (&base_a_$pair and name o6) 2.91 0.01 0.01 
           assign (&base_b_$pair and name n3) 
                  (&base_a_$pair and name n1) 2.95 0.01 0.01
           assign (&base_b_$pair and name o2)
                  (&base_a_$pair and name n2) 2.86 0.01 0.01
           assign (&base_b_$pair and name n3)
                  (&base_a_$pair and name n2) 3.65 0.01 0.01 
           assign (&base_b_$pair and name o2)
                  (&base_a_$pair and name o6) 5.42 0.01 0.01

         end if

       end

     end if

   else

     evaluate ($done=true)

   end if         

   evaluate ($pair=$pair+1)

 end loop noe

 flags include noe end

 set message=$message_old echo=$echo_old end
