1 subroutine read_dist_constr
2 implicit real*8 (a-h,o-z)
4 include 'DIMENSIONS.FREE'
5 include 'COMMON.CONTROL'
7 include 'COMMON.IOUNITS'
8 include 'COMMON.SBRIDGE'
9 integer ifrag_(2,100),ipair_(2,100)
10 double precision wfrag_(100),wpair_(100)
11 character*500 controlcard
12 logical normalize,next
14 double precision xlink(4,0:4) /
16 & 0.0d0,0.0d0,0.0d0,0.0d0, ! default, no xlink potential
17 & 0.00305218d0,9.46638d0,4.68901d0,4.74347d0, ! ZL
18 & 0.00214928d0,12.7517d0,0.00375009d0,6.13477d0, ! ADH
19 & 0.00184547d0,11.2678d0,0.00140292d0,7.00868d0, ! PDH
20 & 0.000161786d0,6.29273d0,4.40993d0,7.13956d0 / ! DSS
21 write (iout,*) "Calling read_dist_constr"
22 c write (iout,*) "nres",nres," nstart_sup",nstart_sup," nsup",nsup
28 call card_concat(controlcard,.true.)
29 next = index(controlcard,"NEXT").gt.0
30 call readi(controlcard,"RESTR_TYPE",restr_type,constr_dist)
31 write (iout,*) "restr_type",restr_type
32 call readi(controlcard,"NFRAG",nfrag_,0)
33 call readi(controlcard,"NPAIR",npair_,0)
34 call readi(controlcard,"NDIST",ndist_,0)
35 call reada(controlcard,'DIST_CUT',dist_cut,5.0d0)
37 & call reada(controlcard,'WBOLTZD',wboltzd,0.591d0)
39 & call reada(controlcard,'SCAL_PEAK',scal_peak,5.0d0)
40 call multreadi(controlcard,"IFRAG",ifrag_(1,1),2*nfrag_,0)
41 call multreadi(controlcard,"IPAIR",ipair_(1,1),2*npair_,0)
42 call multreada(controlcard,"WFRAG",wfrag_(1),nfrag_,0.0d0)
43 call multreada(controlcard,"WPAIR",wpair_(1),npair_,0.0d0)
44 normalize = index(controlcard,"NORMALIZE").gt.0
45 write (iout,*) "WBOLTZD",wboltzd
46 write (iout,*) "SCAL_PEAK",scal_peak
47 c write (iout,*) "NFRAG",nfrag_," NPAIR",npair_," NDIST",ndist_
48 c write (iout,*) "IFRAG"
50 c write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
52 c write (iout,*) "IPAIR"
54 c write (iout,*) i,ipair_(1,i),ipair_(2,i),wpair_(i)
58 read(inp,'(a)') pdbfile
60 & "Distance restraints will be constructed from structure ",pdbfile
61 open(ipdbin,file=pdbfile,status='old',err=11)
67 if (ifrag_(1,i).lt.nstart_sup) ifrag_(1,i)=nstart_sup
68 if (ifrag_(2,i).gt.nstart_sup+nsup-1)
69 & ifrag_(2,i)=nstart_sup+nsup-1
70 c write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
72 if (wfrag_(i).eq.0.0d0) cycle
73 do j=ifrag_(1,i),ifrag_(2,i)-1
75 c write (iout,*) "j",j," k",k
77 if (restr_type.eq.1) then
83 forcon(nhpb)=wfrag_(i)
84 else if (constr_dist.eq.2) then
85 if (ddjk.le.dist_cut) then
91 forcon(nhpb)=wfrag_(i)
93 else if (restr_type.eq.3) then
99 forcon(nhpb)=wfrag_(i)*dexp(-0.5d0*(ddjk/dist_cut)**2)
101 write (iout,'(a,3i5,f8.2,1pe12.2)') "+dist.restr ",
102 & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
107 if (wpair_(i).eq.0.0d0) cycle
115 do j=ifrag_(1,ii),ifrag_(2,ii)
116 do k=ifrag_(1,jj),ifrag_(2,jj)
117 if (restr_type.eq.1) then
123 forcon(nhpb)=wfrag_(i)
124 else if (constr_dist.eq.2) then
125 if (ddjk.le.dist_cut) then
131 forcon(nhpb)=wfrag_(i)
133 else if (restr_type.eq.3) then
139 forcon(nhpb)=wfrag_(i)*dexp(-0.5d0*(ddjk/dist_cut)**2)
141 write (iout,'(a,3i5,f8.2,f10.1)') "+dist.restr ",
142 & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
148 write (iout,*) "Distance restraints as read from input"
150 if (restr_type.eq.12) then
151 read (inp,*) ihpb_peak(nhpb_peak+1),jhpb_peak(nhpb_peak+1),
152 & dhpb_peak(nhpb_peak+1),dhpb1_peak(nhpb_peak+1),
153 & ibecarb_peak(nhpb_peak+1),forcon_peak(nhpb_peak+1),
154 & fordepth_peak(nhpb_peak+1),npeak
155 c write(iout,*) ihpb_peak(nhpb_peak+1),jhpb_peak(nhpb_peak+1),
156 c & dhpb_peak(nhpb_peak+1),dhpb1_peak(nhpb_peak+1),
157 c & ibecarb_peak(nhpb_peak+1),forcon_peak(nhpb_peak+1),
158 c & fordepth_peak(nhpb_peak+1),npeak
159 if (forcon_peak(nhpb_peak+1).le.0.0d0.or.
160 & fordepth_peak(nhpb_peak+1).le.0.0d0)cycle
161 nhpb_peak=nhpb_peak+1
162 irestr_type_peak(nhpb_peak)=12
163 if (ipeak(1,npeak).eq.0) ipeak(1,npeak)=i
165 write (iout,'(a,5i5,2f8.2,2f10.5,i5)') "+dist.restr ",
166 & nhpb_peak,ihpb_peak(nhpb_peak),jhpb_peak(nhpb_peak),
167 & ibecarb_peak(nhpb_peak),npeak,dhpb_peak(nhpb_peak),
168 & dhpb1_peak(nhpb_peak),forcon_peak(nhpb_peak),
169 & fordepth_peak(nhpb_peak),irestr_type_peak(nhpb_peak)
170 if (ibecarb_peak(nhpb_peak).gt.0) then
171 ihpb_peak(nhpb_peak)=ihpb_peak(nhpb_peak)+nres
172 jhpb_peak(nhpb_peak)=jhpb_peak(nhpb_peak)+nres
174 else if (restr_type.eq.11) then
175 read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),dhpb(nhpb+1),
176 & dhpb1(nhpb+1),ibecarb(i),forcon(nhpb+1),fordepth(nhpb+1)
177 c fordepth(nhpb+1)=fordepth(nhpb+1)/forcon(nhpb+1)
178 if (forcon(nhpb+1).le.0.0d0.or.fordepth(nhpb+1).le.0.0d0)cycle
181 write (iout,'(a,4i5,2f8.2,2f10.5,i5)') "+dist.restr ",
182 & nhpb,ihpb(nhpb),jhpb(nhpb),ibecarb(nhpb),dhpb(nhpb),
183 & dhpb1(nhpb),forcon(nhpb),fordepth(nhpb),irestr_type(nhpb)
184 if (ibecarb(nhpb).gt.0) then
185 ihpb(nhpb)=ihpb(nhpb)+nres
186 jhpb(nhpb)=jhpb(nhpb)+nres
188 else if (constr_dist.eq.10) then
189 c Cross-lonk Markov-like potential
190 call card_concat(controlcard,.true.)
191 call readi(controlcard,"ILINK",ihpb(nhpb+1),0)
192 call readi(controlcard,"JLINK",jhpb(nhpb+1),0)
194 if (index(controlcard,"BETA").gt.0) ibecarb(nhpb+1)=1
195 if (ihpb(nhpb+1).eq.0 .or. jhpb(nhpb+1).eq.0) cycle
196 if (index(controlcard,"ZL").gt.0) then
198 else if (index(controlcard,"ADH").gt.0) then
200 else if (index(controlcard,"PDH").gt.0) then
202 else if (index(controlcard,"DSS").gt.0) then
207 call reada(controlcard,"AXLINK",dhpb(nhpb+1),
208 & xlink(1,link_type))
209 call reada(controlcard,"BXLINK",dhpb1(nhpb+1),
210 & xlink(2,link_type))
211 call reada(controlcard,"CXLINK",fordepth(nhpb+1),
212 & xlink(3,link_type))
213 call reada(controlcard,"SIGMA",forcon(nhpb+1),
214 & xlink(4,link_type))
215 call reada(controlcard,"SCORE",xlscore(nhpb+1),1.0d0)
216 c read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),ibecarb(nhpb+1),
217 c & dhpb(nhpb+1),dhpb1(nhpb+1),forcon(nhpb+1),fordepth(nhpb+1)
218 if (forcon(nhpb+1).le.0.0d0 .or.
219 & (dhpb(nhpb+1).eq.0 .and. dhpb1(nhpb+1).eq.0)) cycle
222 if (ibecarb(nhpb).gt.0) then
223 ihpb(nhpb)=ihpb(nhpb)+nres
224 jhpb(nhpb)=jhpb(nhpb)+nres
226 write (iout,'(a,4i5,2f8.2,3f10.5,i5)') "+dist.restr ",
227 & nhpb,ihpb(nhpb),jhpb(nhpb),ibecarb(nhpb),dhpb(nhpb),
228 & dhpb1(nhpb),forcon(nhpb),fordepth(nhpb),xlscore(nhpb),
232 read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),dhpb(nhpb+1),
233 & dhpb1(nhpb+1),ibecarb(nhpb+1),forcon(nhpb+1)
234 if (forcon(nhpb+1).gt.0.0d0) then
236 if (dhpb1(nhpb).eq.0.0d0) then
241 if (ibecarb(nhpb).gt.0) then
242 ihpb(nhpb)=ihpb(nhpb)+nres
243 jhpb(nhpb)=jhpb(nhpb)+nres
245 if (dhpb(nhpb).eq.0.0d0)
246 & dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
248 write (iout,'(a,4i5,f8.2,f10.1)') "+dist.restr ",
249 & nhpb,ihpb(nhpb),jhpb(nhpb),ibecarb(i),dhpb(nhpb),forcon(nhpb)
251 C read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),forcon(nhpb+1)
252 C if (forcon(nhpb+1).gt.0.0d0) then
254 C dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
262 if (irestr_type(i).eq.11.and.fordepth(i).gt.fordepthmax)
263 & fordepthmax=fordepth(i)
266 if (irestr_type(i).eq.11) fordepth(i)=fordepth(i)/fordepthmax
269 if (nhpb.gt.nss) then
270 write (iout,'(/a,i5,a/4a5,2a8,3a10,a5)')
271 & "The following",nhpb-nss,
272 & " distance restraints have been imposed:",
273 & " Nr"," res1"," res2"," beta"," d1"," d2"," k"," V",
276 write (iout,'(4i5,2f8.2,3f10.5,i5)')i-nss,ihpb(i),jhpb(i),
277 & ibecarb(i),dhpb(i),dhpb1(i),forcon(i),fordepth(i),xlscore(i),
285 11 write (iout,*)"read_dist_restr: error reading reference structure"