Merge branch 'lipid' of mmka.chem.univ.gda.pl:unres into lipid
[unres.git] / source / wham / src-NEWSC-NEWCORR / define_pairs.f
diff --git a/source/wham/src-NEWSC-NEWCORR/define_pairs.f b/source/wham/src-NEWSC-NEWCORR/define_pairs.f
new file mode 100644 (file)
index 0000000..00866a8
--- /dev/null
@@ -0,0 +1,120 @@
+      subroutine define_pairs
+      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.FRAG'
+      include 'COMMON.CHAIN'
+      include 'COMMON.HEADER'
+      include 'COMMON.GEO'
+      include 'COMMON.CONTACTS1'
+      include 'COMMON.PEPTCONT'
+      do j=1,nfrag(1)
+        length_frag = 0
+        do k=1,npiece(j,1)
+          length_frag=length_frag+ifrag(2,k,j)-ifrag(1,k,j)+1
+        enddo
+        len_frag(j,1)=length_frag
+        write (iout,*) "Fragment",j," length",len_frag(j,1)
+      enddo
+      nfrag(2)=0
+      do i=1,nfrag(1)
+        do j=i+1,nfrag(1)
+          ind = icant(i,j)
+          if (istruct(i).le.1 .or. istruct(j).le.1) then
+            if (istruct(i).le.1) then
+              ll1=len_frag(i,1)
+            else
+              ll1=len_frag(i,1)/2
+            endif
+            if (istruct(j).le.1) then
+              ll2=len_frag(j,1)
+            else
+              ll2=len_frag(j,1)/2
+            endif
+            len_cut=max0(min0(ll1*2/3,ll2*4/5),3)
+          else
+            if (istruct(i).eq.2 .or. istruct(i).eq.4) then
+              ll1=len_frag(i,1)/2
+            else
+              ll1=len_frag(i,1) 
+            endif
+            if (istruct(j).eq.2 .or. istruct(j).eq.4) then
+              ll2=len_frag(j,1)/2
+            else
+              ll2=len_frag(j,1) 
+            endif
+            len_cut=max0(min0(ll1*4/5,ll2)*4/5,3)
+          endif
+          write (iout,*) "Fragments",i,j," structure",istruct(i),
+     &       istruct(j)," # contacts",
+     &       ncont_frag_ref(ind),nsccont_frag_ref(ind),
+     &       " lengths",len_frag(i,1),len_frag(j,1),
+     &       " ll1",ll1," ll2",ll2," len_cut",len_cut
+          if ((istruct(i).eq.1 .or. istruct(j).eq.1) .and.
+     &      nsccont_frag_ref(ind).ge.len_cut ) then
+            if (istruct(i).eq.1 .and. istruct(j).eq.1) then
+              write (iout,*) "Adding pair of helices",i,j,
+     &        " based on SC contacts"
+            else
+              write (iout,*) "Adding helix+strand/sheet pair",i,j,
+     &        " based on SC contacts"
+            endif
+            nfrag(2)=nfrag(2)+1
+            if (icont_pair.gt.0) then
+              write (iout,*)  "# SC contacts will be used",
+     &        " in comparison."
+              isccont(nfrag(2),2)=1
+            endif
+            if (irms_pair.gt.0) then
+              write (iout,*)  "Fragment RMSD will be used",
+     &        " in comparison."
+              irms(nfrag(2),2)=1
+            endif
+            npiece(nfrag(2),2)=2
+            ipiece(1,nfrag(2),2)=i
+            ipiece(2,nfrag(2),2)=j
+            ielecont(nfrag(2),2)=0
+            n_shift(1,nfrag(2),2)=nshift_pair
+            n_shift(2,nfrag(2),2)=nshift_pair
+            nc_fragm(nfrag(2),2)=ncfrac_pair
+            nc_req_setf(nfrag(2),2)=ncreq_pair
+          else if ((istruct(i).ge.2 .and. istruct(i).le.4)
+     &       .and. (istruct(j).ge.2 .and. istruct(i).le.4)
+     &       .and. ncont_frag_ref(ind).ge.len_cut ) then
+            nfrag(2)=nfrag(2)+1
+            write (iout,*) "Adding pair strands/sheets",i,j,
+     &        " based on pp contacts"
+            if (icont_pair.gt.0) then
+              write (iout,*) "# pp contacts will be used",
+     &        " in comparison."
+              ielecont(nfrag(2),2)=1
+            endif
+            if (irms_pair.gt.0) then
+              write (iout,*)  "Fragment RMSD will be used",
+     &        " in comparison."
+              irms(nfrag(2),2)=1
+            endif
+            npiece(nfrag(2),2)=2
+            ipiece(1,nfrag(2),2)=i
+            ipiece(2,nfrag(2),2)=j
+            ielecont(nfrag(2),2)=1
+            isccont(nfrag(2),2)=0
+            n_shift(1,nfrag(2),2)=nshift_pair
+            n_shift(2,nfrag(2),2)=nshift_pair
+            nc_fragm(nfrag(2),2)=ncfrac_bet
+            nc_req_setf(nfrag(2),2)=ncreq_bet
+          endif
+        enddo
+      enddo
+      write (iout,*) "Pairs found"
+      do i=1,nfrag(2)
+        write (iout,*) ipiece(1,i,2),ipiece(2,i,2)
+      enddo
+      return
+      end