Merge branch 'lipid' of mmka.chem.univ.gda.pl:unres into lipid
[unres.git] / source / wham / src-NEWSC-NEWCORR / proc_cont.f
diff --git a/source/wham/src-NEWSC-NEWCORR/proc_cont.f b/source/wham/src-NEWSC-NEWCORR/proc_cont.f
new file mode 100644 (file)
index 0000000..9269496
--- /dev/null
@@ -0,0 +1,156 @@
+      subroutine proc_cont
+      implicit real*8 (a-h,o-z)
+      include 'DIMENSIONS'
+      include 'DIMENSIONS.ZSCOPT'
+      include 'DIMENSIONS.COMPAR'
+      include 'COMMON.IOUNITS'
+      include 'COMMON.TIME1'
+      include 'COMMON.SBRIDGE'
+      include 'COMMON.CONTROL'
+      include 'COMMON.COMPAR'
+      include 'COMMON.CHAIN'
+      include 'COMMON.HEADER'
+      include 'COMMON.CONTACTS1'
+      include 'COMMON.PEPTCONT'
+      include 'COMMON.GEO'
+      write (iout,*) "proc_cont: nlevel",nlevel
+      if (nlevel.lt.0) then
+        write (iout,*) "call define_fragments"
+        call define_fragments
+      else
+        write (iout,*) "call secondary2"
+        call secondary2(.true.,.false.,ncont_pept_ref,icont_pept_ref,
+     &     isec_ref)
+      endif
+      write (iout,'(80(1h=))')
+      write (iout,*) "Electrostatic contacts"
+      call contacts_between_fragments(.true.,0,ncont_pept_ref,
+     & icont_pept_ref,ncont_frag_ref(1),icont_frag_ref(1,1,1))
+      write (iout,'(80(1h=))')
+      write (iout,*) "Side chain contacts"
+      call contacts_between_fragments(.true.,0,ncont_ref,
+     & icont_ref,nsccont_frag_ref(1),isccont_frag_ref(1,1,1))
+      if (nlevel.lt.0) then
+        do i=1,nfrag(1)
+          ind=icant(i,i)
+          len_cut=1000
+          if (istruct(i).le.1) then
+            len_cut=max0(len_frag(i,1)*4/5,3)
+          else if (istruct(i).eq.2 .or. istruct(i).eq.4) then
+            len_cut=max0(len_frag(i,1)*2/5,3)
+          endif
+          write (iout,*) "i",i," istruct",istruct(i)," ncont_frag",
+     &      ncont_frag_ref(ind)," len_cut",len_cut,
+     &      " icont_single",icont_single," iloc_single",iloc_single
+          iloc(i)=iloc_single
+          if (iloc(i).gt.0) write (iout,*) 
+     &     "Local structure used to compare structure of fragment",i,
+     &     " to native."
+          if (istruct(i).ne.3 .and. istruct(i).ne.0 
+     &        .and. icont_single.gt.0 .and. 
+     &        ncont_frag_ref(ind).ge.len_cut) then
+            write (iout,*) "Electrostatic contacts used to compare",
+     &       " structure of fragment",i," to native."
+            ielecont(i,1)=1
+            isccont(i,1)=0
+          else if (icont_single.gt.0 .and. nsccont_frag_ref(ind)
+     &      .ge.len_cut) then
+            write (iout,*) "Side chain contacts used to compare",
+     &       " structure of fragment",i," to native."
+            isccont(i,1)=1
+            ielecont(i,1)=0
+          else
+            write (iout,*) "Contacts not used to compare",
+     &       " structure of fragment",i," to native."
+            ielecont(i,1)=0
+            isccont(i,1)=0
+            nc_req_setf(i,1)=0
+          endif
+          if (irms_single.gt.0 .or. isccont(i,1).eq.0 
+     &         .and. ielecont(i,1).eq.0) then
+            write (iout,*) "RMSD used to compare",
+     &       " structure of fragment",i," to native."
+            irms(i,1)=1
+          else
+            write (iout,*) "RMSD not used to compare",
+     &       " structure of fragment",i," to native."
+            irms(i,1)=0
+          endif
+        enddo
+      endif
+      if (nlevel.lt.-1) then
+        call define_pairs
+        nlevel = -nlevel
+        if (nlevel.gt.3) nlevel=3
+        if (nlevel.eq.3) then
+          nfrag(3)=1
+          npiece(1,3)=nfrag(1)
+          do i=1,nfrag(1)
+            ipiece(i,1,3)=i
+          enddo
+          ielecont(1,3)=0
+          isccont(1,3)=0
+          irms(1,3)=1
+          n_shift(1,1,3)=0
+          n_shift(2,1,3)=0
+        endif 
+      else if (nlevel.eq.-1) then
+        nlevel=1
+      endif
+      isnfrag(1)=0
+      do i=1,nlevel
+        isnfrag(i+1)=isnfrag(i)+nfrag(i)
+      enddo
+      ndigit=3*nfrag(1)
+      do i=2,nlevel
+        ndigit=ndigit+2*nfrag(i)
+      enddo
+      write (iout,*) "ndigit",ndigit
+      if (.not.binary .and. ndigit.gt.30) then
+        write (iout,*) "Highest class too large; switching to",
+     &    " binary representation."
+        binary=.true.
+      endif
+      write (iout,*) "isnfrag",(isnfrag(i),i=1,nlevel+1)
+      write(iout,*) "rmscut_base_up",rmscut_base_up,
+     & " rmscut_base_low",rmscut_base_low," rmsup_lim",rmsup_lim
+      do i=1,nlevel
+        do j=1,nfrag(i)
+          length_frag = 0
+          if (i.eq.1) then
+            do k=1,npiece(j,i)
+              length_frag=length_frag+ifrag(2,k,j)-ifrag(1,k,j)+1
+            enddo
+          else
+            do k=1,npiece(j,i)
+              length_frag=length_frag+len_frag(ipiece(k,j,i),1)
+            enddo
+          endif
+          len_frag(j,i)=length_frag
+          rmscutfrag(1,j,i)=rmscut_base_up*length_frag
+          rmscutfrag(2,j,i)=rmscut_base_low*length_frag 
+          if (rmscutfrag(1,j,i).lt.rmsup_lim) 
+     &      rmscutfrag(1,j,i)=rmsup_lim
+          if (rmscutfrag(1,j,i).gt.rmsupup_lim) 
+     &      rmscutfrag(1,j,i)=rmsupup_lim
+        enddo
+      enddo
+      write (iout,*) "Level",1," number of fragments:",nfrag(1)
+      do j=1,nfrag(1)
+        write (iout,*) npiece(j,1),(ifrag(1,k,j),ifrag(2,k,j),
+     &    k=1,npiece(j,1)),len_frag(j,1),rmscutfrag(1,j,1),
+     &    rmscutfrag(2,j,1),n_shift(1,j,1),n_shift(2,j,1),
+     &    ang_cut(j)*rad2deg,ang_cut1(j)*rad2deg,frac_min(j),
+     &    nc_fragm(j,1),nc_req_setf(j,1),istruct(j)
+      enddo
+      do i=2,nlevel
+        write (iout,*) "Level",i," number of fragments:",nfrag(i)
+        do j=1,nfrag(i)
+          write (iout,*) npiece(j,i),(ipiece(k,j,i),
+     &      k=1,npiece(j,i)),len_frag(j,i),rmscutfrag(1,j,i),
+     &      rmscutfrag(2,j,i),n_shift(1,j,i),n_shift(2,j,i),
+     &      nc_fragm(j,i),nc_req_setf(j,i) 
+        enddo
+      enddo
+      return
+      end