X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Funres%2Fsrc_MD%2Freadrtns.F;h=b51386c5c5c083bcbbbd37675003613cc902a10c;hb=582d0baf775d8a72f8fd25bac504c98b7f927d86;hp=915e5f23b47aff6ae8ddef197fe07274588713d8;hpb=91c1338cc1805b393e149c3c63a9071f955613ca;p=unres.git diff --git a/source/unres/src_MD/readrtns.F b/source/unres/src_MD/readrtns.F index 915e5f2..b51386c 100644 --- a/source/unres/src_MD/readrtns.F +++ b/source/unres/src_MD/readrtns.F @@ -2581,6 +2581,7 @@ c call flush(iout) call readi(controlcard,"NPAIR",npair_,0) call readi(controlcard,"NDIST",ndist_,0) call reada(controlcard,'DIST_CUT',dist_cut,5.0d0) + call multreadi(controlcard,"IFRAG",ifrag_(1,1),2*nfrag_,0) call multreadi(controlcard,"IPAIR",ipair_(1,1),2*npair_,0) call multreada(controlcard,"WFRAG",wfrag_(1),nfrag_,0.0d0) @@ -2720,6 +2721,7 @@ c------------------------------------------------------------------------------- include 'COMMON.MD' include 'COMMON.GEO' include 'COMMON.INTERACT' + include 'COMMON.NAMES' c c For new homol impl c @@ -2758,7 +2760,7 @@ c Alternative: reading from input call reada(controlcard,"HOMOL_THETA",waga_theta,1.0d0) ! new call reada(controlcard,"HOMOL_SCD",waga_d,1.0d0) ! new call reada(controlcard,'DIST_CUT',dist_cut,5.0d0) ! for diff ways of calc sigma - + call reada(controlcard,'DIST2_CUT',dist2_cut,9999.0d0) call readi(controlcard,"HOMOL_NSET",homol_nset,1) read2sigma=(index(controlcard,'READ2SIGMA').gt.0) start_from_model=(index(controlcard,'START_FROM_MODELS').gt.0) @@ -2793,11 +2795,6 @@ cd call flush(iout) lim_odl=0 lim_dih=0 c -c New -c - lim_theta=0 - lim_xx=0 -c write(iout,*) 'nnt=',nnt,'nct=',nct c do i = nnt,nct @@ -2861,14 +2858,17 @@ c c write(iout,*) "tpl_k_rescore - ",tpl_k_rescore open (ientin,file=tpl_k_rescore,status='old') if (nnt.gt.1) rescore(k,1)=0.0d0 - do irec=nnt,maxdim ! loop for reading res sim + do irec=nnt,nct ! loop for reading res sim if (read2sigma) then read (ientin,*,end=1401) i_tmp,rescore2_tmp,rescore_tmp, & idomain_tmp i_tmp=i_tmp+nnt-1 idomain(k,i_tmp)=idomain_tmp - rescore(k,i_tmp)=0.5d0*(rescore_tmp+0.5d0) - rescore2(k,i_tmp)=0.5d0*(rescore2_tmp+0.5d0) + rescore(k,i_tmp)=rescore_tmp + rescore2(k,i_tmp)=rescore2_tmp + write(iout,'(a7,i5,2f10.5,i5)') "rescore", + & i_tmp,rescore2_tmp,rescore_tmp, + & idomain_tmp else idomain(k,irec)=1 read (ientin,*,end=1401) rescore_tmp @@ -2885,8 +2885,15 @@ c write(iout,*) "rescore(",k,irec,") =",rescore(k,irec) do i = nnt,nct-2 do j=i+2,nct - if (idomain(k,i).eq.idomain(k,j).and.idomain(k,i).ne.0) then - + x12=c(1,i)-c(1,j) + y12=c(2,i)-c(2,j) + z12=c(3,i)-c(3,j) + distal=dsqrt(x12*x12+y12*y12+z12*z12) +c write (iout,*) k,i,j,distal,dist2_cut + + if (idomain(k,i).eq.idomain(k,j).and.idomain(k,i).ne.0 + & .and. distal.le.dist2_cut ) then + ii=ii+1 ii_in_use(ii)=1 l_homo(k,ii)=.true. @@ -2896,13 +2903,6 @@ c write (iout,*) "i",i," j",j," constr_homology", c & constr_homology ires_homo(ii)=i jres_homo(ii)=j -c -c Attempt to replace dist(i,j) by its definition in ... -c - x12=c(1,i)-c(1,j) - y12=c(2,i)-c(2,j) - z12=c(3,i)-c(3,j) - distal=dsqrt(x12*x12+y12*y12+z12*z12) odl(k,ii)=distal if (read2sigma) then sigma_odl(k,ii)=0 @@ -2914,30 +2914,18 @@ c & sigma_odl(k,ii)*dexp(0.5d0*(odl(k,ii)/dist_cut)**2-0.5d0) else 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) + sigma_odl(k,ii)=rescore(k,i)+rescore(k,j) else #ifdef OLDSIGMA - sigma_odl(k,ii)=(rescore(k,i)+rescore(k,j))* ! sigma ~ rescore ~ error + sigma_odl(k,ii)=(rescore(k,i)+rescore(k,j))* & dexp(0.5d0*(odl(k,ii)/dist_cut)**2) #else - sigma_odl(k,ii)=(rescore(k,i)+rescore(k,j))* ! sigma ~ rescore ~ error + sigma_odl(k,ii)=(rescore(k,i)+rescore(k,j))* & 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) -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) -c -c sigma_odl(k,ii)=sigma_odl_temp(i,j,k)* ! new var read from file (?) -c & sigma_odl_temp(i,j,k) ! not inverse because of use of res. similarity + sigma_odl(k,ii)=1.0d0/(sigma_odl(k,ii)*sigma_odl(k,ii)) else ii=ii+1 l_homo(k,ii)=.false. @@ -2980,7 +2968,8 @@ c write (iout,'(i5,10(2f8.2,4x))') i,sigma_dih(k,i) c sigma_dih(k,i)=hmscore(k)*rescore(k,i)*rescore(k,i-1)* c rescore(k,i-2)*rescore(k,i-3) ! right expression ? c Instead of res sim other local measure of b/b str reliability possible - sigma_dih(k,i)=1.0d0/(sigma_dih(k,i)*sigma_dih(k,i)) + if (sigma_dih(k,i).ne.0) + & sigma_dih(k,i)=1.0d0/(sigma_dih(k,i)*sigma_dih(k,i)) c sigma_dih(k,i)=sigma_dih(k,i)*sigma_dih(k,i) enddo lim_dih=nct-nnt-2 @@ -3012,14 +3001,14 @@ c read (ientin,*) sigma_theta(k,i) ! 1st variant sigma_theta(k,i)=(rescore(k,i)+rescore(k,i-1)+ & rescore(k,i-2))/3.0 c if (read2sigma) sigma_theta(k,i)=sigma_theta(k,i)/3.0 - sigma_theta(k,i)=1.0d0/(sigma_theta(k,i)*sigma_theta(k,i)) + if (sigma_theta(k,i).ne.0) + & sigma_theta(k,i)=1.0d0/(sigma_theta(k,i)*sigma_theta(k,i)) c sigma_theta(k,i)=hmscore(k)*rescore(k,i)*rescore(k,i-1)* c rescore(k,i-2) ! right expression ? c sigma_theta(k,i)=sigma_theta(k,i)*sigma_theta(k,i) enddo endif - lim_theta=nct-nnt-1 if (waga_d.gt.0.0d0) then c open (ientin,file=tpl_k_sigma_d,status='old') @@ -3046,14 +3035,14 @@ c write (iout,*) "yytpl(",k,i,") =",yytpl(k,i) c write (iout,*) "zztpl(",k,i,") =",zztpl(k,i) c write(iout,*) "rescore(",k,i,") =",rescore(k,i) sigma_d(k,i)=rescore(k,i) ! right expression ? - sigma_d(k,i)=1.0d0/(sigma_d(k,i)*sigma_d(k,i)) + if (sigma_d(k,i).ne.0) + & sigma_d(k,i)=1.0d0/(sigma_d(k,i)*sigma_d(k,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? enddo - lim_xx=nct-nnt+1 endif enddo c @@ -3096,24 +3085,26 @@ cd write(iout,*) "waga_theta",waga_theta,"waga_d",waga_d if(me.eq.king .or. .not. out1file .and. fg_rank.eq.0) then write (iout,*) "Distance restraints from templates" do ii=1,lim_odl - write(iout,'(3i5,10(2f8.2,1x,l1,4x))') + write(iout,'(3i5,100(2f8.2,1x,l1,4x))') & ii,ires_homo(ii),jres_homo(ii), & (odl(ki,ii),1.0d0/dsqrt(sigma_odl(ki,ii)),l_homo(ki,ii), & ki=1,constr_homology) enddo write (iout,*) "Dihedral angle restraints from templates" - do i=nnt+3,lim_dih - write (iout,'(i5,10(2f8.2,4x))') i,(rad2deg*dih(ki,i), + do i=nnt+3,nct + write (iout,'(i5,a4,100(2f8.2,4x))') i,restyp(itype(i)), + & (rad2deg*dih(ki,i), & rad2deg/dsqrt(sigma_dih(ki,i)),ki=1,constr_homology) enddo write (iout,*) "Virtual-bond angle restraints from templates" - do i=nnt+2,lim_theta - write (iout,'(i5,10(2f8.2,4x))') i,(rad2deg*thetatpl(ki,i), + do i=nnt+2,nct + write (iout,'(i5,a4,100(2f8.2,4x))') i,restyp(itype(i)), + & (rad2deg*thetatpl(ki,i), & rad2deg/dsqrt(sigma_theta(ki,i)),ki=1,constr_homology) enddo write (iout,*) "SC restraints from templates" - do i=nnt,lim_xx - write(iout,'(i5,10(4f8.2,4x))') i, + do i=nnt,nct + write(iout,'(i5,100(4f8.2,4x))') i, & (xxtpl(ki,i),yytpl(ki,i),zztpl(ki,i), & 1.0d0/dsqrt(sigma_d(ki,i)),ki=1,constr_homology) enddo