read2sigma in wham and cluster_wham
[unres.git] / source / wham / src-M / define_pairs.f
1       subroutine define_pairs
2       implicit real*8 (a-h,o-z)
3       include 'DIMENSIONS'
4       include 'DIMENSIONS.ZSCOPT'
5       include 'DIMENSIONS.COMPAR'
6       include 'COMMON.IOUNITS'
7       include 'COMMON.TIME1'
8       include 'COMMON.SBRIDGE'
9       include 'COMMON.CONTROL'
10       include 'COMMON.COMPAR'
11       include 'COMMON.FRAG'
12       include 'COMMON.CHAIN'
13       include 'COMMON.HEADER'
14       include 'COMMON.GEO'
15       include 'COMMON.CONTACTS1'
16       include 'COMMON.PEPTCONT'
17       do j=1,nfrag(1)
18         length_frag = 0
19         do k=1,npiece(j,1)
20           length_frag=length_frag+ifrag(2,k,j)-ifrag(1,k,j)+1
21         enddo
22         len_frag(j,1)=length_frag
23         write (iout,*) "Fragment",j," length",len_frag(j,1)
24       enddo
25       nfrag(2)=0
26       do i=1,nfrag(1)
27         do j=i+1,nfrag(1)
28           ind = icant(i,j)
29           if (istruct(i).le.1 .or. istruct(j).le.1) then
30             if (istruct(i).le.1) then
31               ll1=len_frag(i,1)
32             else
33               ll1=len_frag(i,1)/2
34             endif
35             if (istruct(j).le.1) then
36               ll2=len_frag(j,1)
37             else
38               ll2=len_frag(j,1)/2
39             endif
40             len_cut=max0(min0(ll1*2/3,ll2*4/5),3)
41           else
42             if (istruct(i).eq.2 .or. istruct(i).eq.4) then
43               ll1=len_frag(i,1)/2
44             else
45               ll1=len_frag(i,1) 
46             endif
47             if (istruct(j).eq.2 .or. istruct(j).eq.4) then
48               ll2=len_frag(j,1)/2
49             else
50               ll2=len_frag(j,1) 
51             endif
52             len_cut=max0(min0(ll1*4/5,ll2)*4/5,3)
53           endif
54           write (iout,*) "Fragments",i,j," structure",istruct(i),
55      &       istruct(j)," # contacts",
56      &       ncont_frag_ref(ind),nsccont_frag_ref(ind),
57      &       " lengths",len_frag(i,1),len_frag(j,1),
58      &       " ll1",ll1," ll2",ll2," len_cut",len_cut
59           if ((istruct(i).eq.1 .or. istruct(j).eq.1) .and.
60      &      nsccont_frag_ref(ind).ge.len_cut ) then
61             if (istruct(i).eq.1 .and. istruct(j).eq.1) then
62               write (iout,*) "Adding pair of helices",i,j,
63      &        " based on SC contacts"
64             else
65               write (iout,*) "Adding helix+strand/sheet pair",i,j,
66      &        " based on SC contacts"
67             endif
68             nfrag(2)=nfrag(2)+1
69             if (icont_pair.gt.0) then
70               write (iout,*)  "# SC contacts will be used",
71      &        " in comparison."
72               isccont(nfrag(2),2)=1
73             endif
74             if (irms_pair.gt.0) then
75               write (iout,*)  "Fragment RMSD will be used",
76      &        " in comparison."
77               irms(nfrag(2),2)=1
78             endif
79             npiece(nfrag(2),2)=2
80             ipiece(1,nfrag(2),2)=i
81             ipiece(2,nfrag(2),2)=j
82             ielecont(nfrag(2),2)=0
83             n_shift(1,nfrag(2),2)=nshift_pair
84             n_shift(2,nfrag(2),2)=nshift_pair
85             nc_fragm(nfrag(2),2)=ncfrac_pair
86             nc_req_setf(nfrag(2),2)=ncreq_pair
87           else if ((istruct(i).ge.2 .and. istruct(i).le.4)
88      &       .and. (istruct(j).ge.2 .and. istruct(i).le.4)
89      &       .and. ncont_frag_ref(ind).ge.len_cut ) then
90             nfrag(2)=nfrag(2)+1
91             write (iout,*) "Adding pair strands/sheets",i,j,
92      &        " based on pp contacts"
93             if (icont_pair.gt.0) then
94               write (iout,*) "# pp contacts will be used",
95      &        " in comparison."
96               ielecont(nfrag(2),2)=1
97             endif
98             if (irms_pair.gt.0) then
99               write (iout,*)  "Fragment RMSD will be used",
100      &        " in comparison."
101               irms(nfrag(2),2)=1
102             endif
103             npiece(nfrag(2),2)=2
104             ipiece(1,nfrag(2),2)=i
105             ipiece(2,nfrag(2),2)=j
106             ielecont(nfrag(2),2)=1
107             isccont(nfrag(2),2)=0
108             n_shift(1,nfrag(2),2)=nshift_pair
109             n_shift(2,nfrag(2),2)=nshift_pair
110             nc_fragm(nfrag(2),2)=ncfrac_bet
111             nc_req_setf(nfrag(2),2)=ncreq_bet
112           endif
113         enddo
114       enddo
115       write (iout,*) "Pairs found"
116       do i=1,nfrag(2)
117         write (iout,*) ipiece(1,i,2),ipiece(2,i,2)
118       enddo
119       return
120       end