character*500 controlcard
integer ki, i, j, k, l, ii_in_use(maxdim)
logical lprn /.true./
+ integer ilen
+ external ilen
c
c FP - Nov. 2014 Temporary specifications for new vars
c
- double precision rescore_tmp,x12,y12,z12,rescore2_tmp
+ double precision rescore_tmp,x12,y12,z12
double precision, dimension (max_template,maxres) :: rescore
- double precision, dimension (max_template,maxres) :: rescore2
character*24 pdbfile,tpl_k_rescore
c -----------------------------------------------------------------
c Reading multiple PDB ref structures and calculation of retraints
call reada(controlcard,'DIST_CUT',dist_cut,5.0d0) ! for diff ways of calc sigma
call readi(controlcard,"HOMOL_NSET",homol_nset,1)
- read2sigma=(index(controlcard,'READ2SIGMA').gt.0)
if (homol_nset.gt.1)then
call card_concat(controlcard)
read(controlcard,*) (waga_homology(i),i=1,homol_nset)
lim_theta=0
lim_xx=0
c
+c Reading HM global scores (prob not required)
+c
+c open (4,file="HMscore")
+c do k=1,constr_homology
+c read (4,*,end=521) hmscore_tmp
+c hmscore(k)=hmscore_tmp ! Another transformation can be used
+c write(*,*) "Model", k, ":", hmscore(k)
+c enddo
+c521 continue
+
+c write(iout,*) "waga_theta",waga_theta,"waga_d",waga_d
ii=0
do i = nnt,nct-2
do j=i+2,nct
do k=1,constr_homology
+ read(inp,'(a)') pdbfile
+c Next stament causes error upon compilation (?)
+c if(me.eq.king.or. .not. out1file)
+c write (iout,'(2a)') 'PDB data will be read from file ',
+c & pdbfile(:ilen(pdbfile))
+ write (iout,'(a,5x,a)') 'HOMOL: Opening PDB file',
+ & pdbfile(:ilen(pdbfile))
+ open(ipdbin,file=pdbfile,status='old',err=33)
+ goto 34
+ 33 write (iout,'(a,5x,a)') 'Error opening PDB file',
+ & pdbfile(:ilen(pdbfile))
+ stop
+ 34 continue
+c print *,'Begin reading pdb data'
c
c Files containing res sim or local scores (former containing sigmas)
c
write(kic2,'(bz,i2.2)') k
tpl_k_rescore="template"//kic2//".sco"
+c tpl_k_sigma_odl="template"//kic2//".sigma_odl"
+c tpl_k_sigma_dih="template"//kic2//".sigma_dih"
+c tpl_k_sigma_theta="template"//kic2//".sigma_theta"
+c tpl_k_sigma_d="template"//kic2//".sigma_d"
+ unres_pdb=.false.
+ call readpdb
+c
+c Distance restraints
+c
+c ... --> odl(k,ii)
+C Copy the coordinates from reference coordinates (?)
+ do i=1,2*nres
+ do j=1,3
+ c(j,i)=cref(j,i)
+c write (iout,*) "c(",j,i,") =",c(j,i)
+ enddo
+ enddo
c
c From read_dist_constr (commented out 25/11/2014 <-> res sim)
c
rescore(k,irec)=0.0d0
goto 1301
endif
- if (read2sigma) then
- read (ientin,*,end=1401) rescore2_tmp,rescore_tmp,
- & idomain(k,irec)
- rescore(k,irec)=0.5d0*(rescore_tmp+0.5d0)
- rescore2(k,irec)=0.5d0*(rescore2_tmp+0.5d0)
- else
- idomain(k,irec)=1
- read (ientin,*,end=1401) rescore_tmp
-
+ read (ientin,*,end=1401) rescore_tmp
c rescore(k,irec)=rescore_tmp+1.0d0 ! to avoid 0 values
- rescore(k,irec)=0.5d0*(rescore_tmp+0.5d0) ! alt transf to reduce scores
+ rescore(k,irec)=0.5d0*(rescore_tmp+0.5d0) ! alt transf to reduce scores
c write(iout,*) "rescore(",k,irec,") =",rescore(k,irec)
- endif
1301 continue
enddo
1401 continue
close (ientin)
-
- read(inp,'(a)') pdbfile
-c Next stament causes error upon compilation (?)
-c if(me.eq.king.or. .not. out1file)
-c write (iout,'(2a)') 'PDB data will be read from file ',
-c & pdbfile(:ilen(pdbfile))
- open(ipdbin,file=pdbfile,status='old',err=33)
- goto 34
- 33 write (iout,'(a)') 'Error opening PDB file.'
- stop
- 34 continue
-c print *,'Begin reading pdb data'
-
- unres_pdb=.false.
- if (read2sigma) then
- call readpdb_template(k)
- else
- call readpdb
- endif
-c
-C Copy the coordinates from reference coordinates (?)
- do i=1,2*nres
- do j=1,3
- c(j,i)=cref(j,i)
-c write (iout,*) "c(",j,i,") =",c(j,i)
- enddo
- enddo
-
-
-
+c open (ientin,file=tpl_k_sigma_odl,status='old')
+c do irec=1,maxdim ! loop for reading sigma_odl
+c read (ientin,*,end=1401) i, j,
+c & sigma_odl_temp(i+nnt-1,j+nnt-1,k) ! new variable (?)
+c sigma_odl_temp(j+nnt-1,i+nnt-1,k)= ! which purpose?
+c & sigma_odl_temp(i+nnt-1,j+nnt-1,k)
+c enddo
+c 1401 continue
+c close (ientin)
if (waga_dist.ne.0.0d0) then
ii=0
do i = nnt,nct-2 ! right? without parallel.
c sigma_odl(k,ii)=hmscore(k)*rescore(k,i)*rescore(k,j)
c Other functional forms possible depending on odl(k,ii), eg.
c
- if (read2sigma) then
- sigma_odl(k,ii)=0
- do ik=i,j
- sigma_odl(k,ii)=sigma_odl(k,ii)+rescore2(k,ik)
- enddo
- sigma_odl(k,ii)=sigma_odl(k,ii)/(j-i+1)
- if (odl(k,ii).gt.dist_cut) sigma_odl(k,ii) =
- & sigma_odl(k,ii)*dexp(0.5d0*(odl(k,ii)/dist_cut)**2-0.5d0)
- else
- if (odl(k,ii).le.dist_cut) then
+ if (odl(k,ii).le.dist_cut) then
sigma_odl(k,ii)=rescore(k,i)+rescore(k,j) ! other exprs possible
c sigma_odl(k,ii)=hmscore(k)*rescore(k,i)*rescore(k,j)
else
sigma_odl(k,ii)=(rescore(k,i)+rescore(k,j))* ! sigma ~ rescore ~ error
& dexp(0.5d0*(odl(k,ii)/dist_cut)**2-0.5d0)
#endif
- endif
- endif
+
c Following expr replaced by a positive exp argument
c sigma_odl(k,ii)=(rescore(k,i)+rescore(k,j))*
c & dexp(-0.5d0*(odl(k,ii)/dist_cut)**2)
c sigma_odl(k,ii)=hmscore(k)*rescore(k,i)*rescore(k,j)*
c & dexp(-0.5d0*(odl(k,ii)/dist_cut)**2)
-
+ endif
c
sigma_odl(k,ii)=1.0d0/(sigma_odl(k,ii)*sigma_odl(k,ii)) ! rescore ~ error
c sigma_odl(k,ii)=sigma_odl(k,ii)*sigma_odl(k,ii)
l_homo(k,ii)=.false.
endif
enddo
+c read (ientin,*) sigma_odl(k,ii) ! 1st variant
enddo
+c lim_odl=ii
+c if (constr_homology.gt.0) call homology_partition
lim_odl=ii
endif
c
c do i=1,nres ! alternative for bounds acc to readpdb?
c do i=1,nres-3 ! alternative for bounds as used to set initial values in orig. read_constr_homology
c do i=idihconstr_start_homo,idihconstr_end_homo ! with FG parallel.
- if (idomain(k,i).eq.0) then
- sigma_dih(k,i)=0.0
- cycle
- endif
dih(k,i)=phiref(i) ! right?
c read (ientin,*) sigma_dih(k,i) ! original variant
c write (iout,*) "dih(",k,i,") =",dih(k,i)
do i = nnt+2,nct ! right? without parallel.
c do i = i=1,nres ! alternative for bounds acc to readpdb?
c do i=ithet_start,ithet_end ! with FG parallel.
- if (idomain(k,i).eq.0) then
- sigma_theta(k,i)=0.0
- cycle
- endif
thetatpl(k,i)=thetaref(i)
c write (iout,*) "thetatpl(",k,i,") =",thetatpl(k,i)
c write(iout,*) "rescore(",k,i,") =",rescore(k,i),
c & sigma_d(k,i+nnt-1)
c enddo
c1404 continue
-
+ close (ientin)
do i = nnt,nct ! right? without parallel.
c do i=2,nres-1 ! alternative for bounds acc to readpdb?
c do i=loc_start,loc_end ! with FG parallel.
- if (itype(i).eq.10) cycle
- if (idomain(k,i).eq.0 ) then
- sigma_d(k,i)=0.0
- cycle
- endif
+ if (itype(i).eq.10) goto 1 ! right?
xxtpl(k,i)=xxref(i)
yytpl(k,i)=yyref(i)
zztpl(k,i)=zzref(i)
c sigma_d(k,i)=hmscore(k)*rescore(k,i) ! right expression ?
c sigma_d(k,i)=sigma_d(k,i)*sigma_d(k,i)
c read (ientin,*) sigma_d(k,i) ! 1st variant
+ if (i-nnt+1.gt.lim_xx) lim_xx=i-nnt+1 ! right?
+ 1 continue
enddo
lim_xx=nct-nnt+1
endif
+ close(ientin)
enddo
c
c remove distance restraints not used in any model from the list
projects / unres.git / blobdiff