5 implicit real*8 (a-h,o-z)
7 include 'DIMENSIONS.ZSCOPT'
8 include 'COMMON.IOUNITS'
11 include 'COMMON.INTERACT'
12 include 'COMMON.LOCAL'
13 include 'COMMON.NAMES'
14 include 'COMMON.CHAIN'
15 include 'COMMON.FFIELD'
16 include 'COMMON.SBRIDGE'
17 include 'COMMON.TORCNSTR'
18 include 'COMMON.CONTROL'
19 character*4 sequence(maxres)
21 double precision x(maxvar)
22 character*320 controlcard,ucase
23 dimension itype_pdb(maxres)
25 call card_concat(controlcard,.true.)
26 call reada(controlcard,'SCAL14',scal14,0.4d0)
27 call reada(controlcard,'SCALSCP',scalscp,1.0d0)
28 call reada(controlcard,'CUTOFF',cutoff_corr,7.0d0)
29 call reada(controlcard,'DELT_CORR',delt_corr,0.5d0)
30 r0_corr=cutoff_corr-delt_corr
31 call readi(controlcard,"NRES",nres,0)
32 iscode=index(controlcard,"ONE_LETTER")
34 write (iout,*) "Error: no residues in molecule"
37 if (nres.gt.maxres) then
38 write (iout,*) "Error: too many residues",nres,maxres
40 write(iout,*) 'nres=',nres
41 C Read sequence of the protein
43 read (inp,'(80a1)') (sequence(i)(1:1),i=1,nres)
45 read (inp,'(20(1x,a3))') (sequence(i),i=1,nres)
47 C Convert sequence to numeric code
49 itype(i)=rescode(i,sequence(i),iscode)
51 write (iout,*) "Numeric code:"
52 write (iout,'(20i4)') (itype(i),i=1,nres)
55 if (itype(i).eq.21 .or. itype(i+1).eq.21) then
57 if (itype(i).eq.21) then
61 else if (itype(i+1).ne.20) then
63 else if (itype(i).ne.20) then
72 if (with_dihed_constr) then
74 read (inp,*) ndih_constr
75 if (ndih_constr.gt.0) then
77 write (iout,*) 'FTORS',ftors
78 read (inp,*) (idih_constr(i),phi0(i),drange(i),i=1,ndih_constr)
80 & 'There are',ndih_constr,' constraints on phi angles.'
82 write (iout,'(i5,2f8.3)') idih_constr(i),phi0(i),drange(i)
85 phi0(i)=deg2rad*phi0(i)
86 drange(i)=deg2rad*drange(i)
94 if (itype(1).eq.21) nnt=2
95 if (itype(nres).eq.21) nct=nct-1
96 write(iout,*) 'NNT=',NNT,' NCT=',NCT
97 c Read distance restraints
98 if (constr_dist.gt.0) then
99 if (refstr) call read_ref_structure(*11)
100 call read_dist_constr
104 if (constr_homology.gt.0) then
105 call read_constr_homology
112 write (iout,'(/a,i3,a)') 'The chain contains',ns,
113 & ' disulfide-bridging cysteines.'
114 write (iout,'(20i4)') (iss(i),i=1,ns)
115 write (iout,'(/a/)') 'Pre-formed links are:'
121 write (iout,'(2a,i3,3a,i3,a,3f10.3)')
122 & restyp(it1),'(',i1,') -- ',restyp(it2),'(',i2,')',
123 & dhpb(i),ebr,forcon(i)
128 11 stop "Error reading reference structure"
130 c-----------------------------------------------------------------------------
131 logical function seq_comp(itypea,itypeb,length)
133 integer length,itypea(length),itypeb(length)
136 if (itypea(i).ne.itypeb(i)) then
144 c-----------------------------------------------------------------------------
145 subroutine read_bridge
146 C Read information about disulfide bridges.
147 implicit real*8 (a-h,o-z)
149 include 'DIMENSIONS.ZSCOPT'
150 include 'COMMON.IOUNITS'
153 include 'COMMON.INTERACT'
154 include 'COMMON.NAMES'
155 include 'COMMON.CHAIN'
156 include 'COMMON.FFIELD'
157 include 'COMMON.SBRIDGE'
158 C Read bridging residues.
159 read (inp,*) ns,(iss(i),i=1,ns)
161 write (iout,*) 'ns=',ns,' iss:',(iss(i),i=1,ns)
162 C Check whether the specified bridging residues are cystines.
164 if (itype(iss(i)).ne.1) then
165 write (iout,'(2a,i3,a)')
166 & 'Do you REALLY think that the residue ',restyp(iss(i)),i,
167 & ' can form a disulfide bridge?!!!'
168 write (*,'(2a,i3,a)')
169 & 'Do you REALLY think that the residue ',restyp(iss(i)),i,
170 & ' can form a disulfide bridge?!!!'
174 C Read preformed bridges.
176 read (inp,*) nss,(ihpb(i),jhpb(i),i=1,nss)
177 write (iout,*) 'nss=',nss,' ihpb,jhpb: ',(ihpb(i),jhpb(i),i=1,nss)
180 C Check if the residues involved in bridges are in the specified list of
184 if (ihpb(i).eq.ihpb(j).or.ihpb(i).eq.jhpb(j)
185 & .or.jhpb(i).eq.ihpb(j).or.jhpb(i).eq.jhpb(j)) then
186 write (iout,'(a,i3,a)') 'Disulfide pair',i,
187 & ' contains residues present in other pairs.'
188 write (*,'(a,i3,a)') 'Disulfide pair',i,
189 & ' contains residues present in other pairs.'
194 if (ihpb(i).eq.iss(j)) goto 10
196 write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
199 if (jhpb(i).eq.iss(j)) goto 20
201 write (iout,'(a,i3,a)') 'Pair',i,' contains unknown cystine.'
212 if (ns.gt.0.and.dyn_ss) then
213 C /06/28/2013 Adasko:ns is number of Cysteins bonded also called half of
216 C /06/28/2013 Adasko: nss number of full SS bonds
219 forcon(i-nss)=forcon(i)
226 dyn_ss_mask(iss(i))=.true.
227 C /06/28/2013 Adasko: dyn_ss_mask which Cysteins can form disulfidebond
228 c write(iout,*) i,iss(i),dyn_ss_mask(iss(i)),"ATU"
233 c------------------------------------------------------------------------------
234 subroutine read_angles(kanal,iscor,energ,iprot,*)
235 implicit real*8 (a-h,o-z)
237 include 'DIMENSIONS.ZSCOPT'
238 include 'COMMON.INTERACT'
239 include 'COMMON.SBRIDGE'
242 include 'COMMON.CHAIN'
243 include 'COMMON.IOUNITS'
245 read(kanal,'(a80)',end=10,err=10) lineh
246 read(lineh(:5),*,err=8) ic
247 read(lineh(6:),*,err=8) energ
250 print *,'error, assuming e=1d10',lineh
254 read(lineh(18:),*,end=10,err=10) nss
256 read (lineh(20:),*,end=10,err=10)
257 & (IHPB(I),JHPB(I),I=1,NSS),iscor
259 read (lineh(20:),*,end=10,err=10) (IHPB(I),JHPB(I),I=1,8)
260 read (kanal,*,end=10,err=10) (IHPB(I),JHPB(I),
263 c print *,"energy",energ," iscor",iscor
264 read (kanal,*,err=10,end=10) (theta(i),i=3,nres)
265 read (kanal,*,err=10,end=10) (phi(i),i=4,nres)
266 read (kanal,*,err=10,end=10) (alph(i),i=2,nres-1)
267 read (kanal,*,err=10,end=10) (omeg(i),i=2,nres-1)
269 theta(i)=deg2rad*theta(i)
270 phi(i)=deg2rad*phi(i)
271 alph(i)=deg2rad*alph(i)
272 omeg(i)=deg2rad*omeg(i)
277 c-------------------------------------------------------------------------------
278 subroutine read_dist_constr
279 implicit real*8 (a-h,o-z)
281 include 'COMMON.CONTROL'
282 include 'COMMON.CHAIN'
283 include 'COMMON.IOUNITS'
284 include 'COMMON.SBRIDGE'
285 integer ifrag_(2,100),ipair_(2,100)
286 double precision wfrag_(100),wpair_(100)
287 character*500 controlcard
288 c write (iout,*) "Calling read_dist_constr"
289 c write (iout,*) "nres",nres," nstart_sup",nstart_sup," nsup",nsup
291 call card_concat(controlcard,.true.)
292 call readi(controlcard,"NFRAG",nfrag_,0)
293 call readi(controlcard,"NPAIR",npair_,0)
294 call readi(controlcard,"NDIST",ndist_,0)
295 call reada(controlcard,'DIST_CUT',dist_cut,5.0d0)
296 call multreadi(controlcard,"IFRAG",ifrag_(1,1),2*nfrag_,0)
297 call multreadi(controlcard,"IPAIR",ipair_(1,1),2*npair_,0)
298 call multreada(controlcard,"WFRAG",wfrag_(1),nfrag_,0.0d0)
299 call multreada(controlcard,"WPAIR",wpair_(1),npair_,0.0d0)
300 write (iout,*) "NFRAG",nfrag_," NPAIR",npair_," NDIST",ndist_
301 write (iout,*) "IFRAG"
303 write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
305 write (iout,*) "IPAIR"
307 write (iout,*) i,ipair_(1,i),ipair_(2,i),wpair_(i)
310 if (.not.refstr .and. nfrag_.gt.0) then
312 & "ERROR: no reference structure to compute distance restraints"
314 & "Restraints must be specified explicitly (NDIST=number)"
317 if (nfrag_.lt.2 .and. npair_.gt.0) then
318 write (iout,*) "ERROR: Less than 2 fragments specified",
319 & " but distance restraints between pairs requested"
324 if (ifrag_(1,i).lt.nstart_sup) ifrag_(1,i)=nstart_sup
325 if (ifrag_(2,i).gt.nstart_sup+nsup-1)
326 & ifrag_(2,i)=nstart_sup+nsup-1
327 c write (iout,*) i,ifrag_(1,i),ifrag_(2,i),wfrag_(i)
329 if (wfrag_(i).gt.0.0d0) then
330 do j=ifrag_(1,i),ifrag_(2,i)-1
332 write (iout,*) "j",j," k",k
334 if (constr_dist.eq.1) then
339 forcon(nhpb)=wfrag_(i)
340 else if (constr_dist.eq.2) then
341 if (ddjk.le.dist_cut) then
346 forcon(nhpb)=wfrag_(i)
353 forcon(nhpb)=wfrag_(i)*dexp(-0.5d0*(ddjk/dist_cut)**2)
355 write (iout,'(a,3i5,f8.2,1pe12.2)') "+dist.constr ",
356 & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
362 if (wpair_(i).gt.0.0d0) then
370 do j=ifrag_(1,ii),ifrag_(2,ii)
371 do k=ifrag_(1,jj),ifrag_(2,jj)
375 forcon(nhpb)=wpair_(i)
377 write (iout,'(a,3i5,f8.2,f10.1)') "+dist.constr ",
378 & nhpb,ihpb(nhpb),jhpb(nhpb),dhpb(nhpb),forcon(nhpb)
384 read (inp,*) ihpb(nhpb+1),jhpb(nhpb+1),dhpb(i),dhpb1(i),
385 & ibecarb(i),forcon(nhpb+1)
386 if (forcon(nhpb+1).gt.0.0d0) then
388 if (ibecarb(i).gt.0) then
392 if (dhpb(nhpb).eq.0.0d0)
393 & dhpb(nhpb)=dist(ihpb(nhpb),jhpb(nhpb))
397 write (iout,'(a,3i5,2f8.2,i2,f10.1)') "+dist.constr ",
398 & i,ihpb(i),jhpb(i),dhpb(i),dhpb1(i),ibecarb(i),forcon(i)
406 c====-------------------------------------------------------------------
408 subroutine read_constr_homology
414 c include 'COMMON.SETUP'
415 include 'COMMON.CONTROL'
416 include 'COMMON.CHAIN'
417 include 'COMMON.IOUNITS'
418 c include 'COMMON.MD'
421 character*24 model_ki_dist, model_ki_angle
422 character*500 controlcard
423 integer ki, i, j, k, l
425 c write(iout,*)"TEST1", waga_dist, waga_angle
426 c call card_concat(controlcard)
427 call card_concat(controlcard,.true.)
428 c write(iout,*)"TEST1a", waga_dist, waga_angle
429 call reada(controlcard,"HOMOL_DIST",waga_dist,1.0)
430 c write(iout,*)"TEST1b", waga_dist, waga_angle
431 call reada(controlcard,"HOMOL_ANGLE",waga_angle,1.0)
433 c write(iout,*)"TEST2", waga_dist, waga_angle
435 c write(iout,*) "TEST_constr_homology_czytanie", constr_homology
436 do ki=1,constr_homology
437 if (constr_homology.ge.1) then
439 write(kic2,'(i2)') ki
440 c write(iout,*) "TEST KICA, HOMOL", kic2
441 if (ki.le.9) kic2="0"//kic2(2:2)
442 c write(iout,*) "TEST KICA2, HOMOL", kic2
443 model_ki_dist="model"//kic2//".dist"
444 model_ki_angle="model"//kic2//".angle"
445 c write(iout,*) model_ki_dist, model_ki_angle
446 open (1400+ki,file=model_ki_dist,status='old')
447 open (1401+ki,file=model_ki_angle,status='old')
449 do irec=1,99999 !petla do czytania wiezow na odleglosc
450 read (1400+ki,*,end=1401) i, j, odl(i,j,ki),sigma_odl(i,j,ki)
454 c write(iout,*) "TEST_CZYTANIA__WIEZOW_ODL", i, j, odl(i,j,ki)
455 do irec=1,99999 !petla do czytania wiezow na katach torsyjnych
456 read (1401+ki,*,end=1402) i, j, k,l,dih(i,ki),sigma_dih(i,ki)
458 c dih(i,ki)=dih(i,ki)
464 c write(iout,*) "TEST CZYTANIA1",odl(1,2,1),odl(1,3,1),odl(1,4,1)
465 c write(iout,*) "TEST CZYTANIA2",dih(1,1),dih(2,1),dih(3,1)