From: Pawel Krupa Date: Thu, 10 Mar 2016 22:32:05 +0000 (+0100) Subject: Calculation of folded fractions (at given RMSD cutoff) in WHAM analysis (7th column... X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=commitdiff_plain;h=69f9281377be5253df03c0af9087353950e07aaf;p=unres.git Calculation of folded fractions (at given RMSD cutoff) in WHAM analysis (7th column in the thermal files). --- diff --git a/source/wham/src-M/COMMON.ENERGIES b/source/wham/src-M/COMMON.ENERGIES index 2d40a95..1c82395 100644 --- a/source/wham/src-M/COMMON.ENERGIES +++ b/source/wham/src-M/COMMON.ENERGIES @@ -1,4 +1,4 @@ double precision potE(MaxStr_Proc,Max_Parm),entfac(MaxStr_Proc), - & q(MaxQ+2,MaxStr_Proc),enetb(max_ene,MaxStr_Proc,Max_Parm) + & q(MaxQ+3,MaxStr_Proc),enetb(max_ene,MaxStr_Proc,Max_Parm) integer einicheck common /energies/ potE,entfac,q,enetb,einicheck diff --git a/source/wham/src-M/COMMON.PROTFILES b/source/wham/src-M/COMMON.PROTFILES index 3287326..f5a1b80 100644 --- a/source/wham/src-M/COMMON.PROTFILES +++ b/source/wham/src-M/COMMON.PROTFILES @@ -5,6 +5,7 @@ & nfile_cx(MaxR,MaxT_h,Max_Parm), & rec_start(MaxR,MaxT_h,Max_Parm), & rec_end(MaxR,MaxT_h,Max_Parm),lenrec,lenrec1,lenrec2 - common /protfil/ protfiles,bprotfiles, + real*8 frac_cutoff + common /protfil/ protfiles,bprotfiles, frac_cutoff, & nfile_bin,nfile_asc,nfile_cx,rec_start,rec_end,lenrec,lenrec1, & lenrec2 diff --git a/source/wham/src-M/DIMENSIONS.FREE b/source/wham/src-M/DIMENSIONS.FREE index 691d9b2..304e83d 100644 --- a/source/wham/src-M/DIMENSIONS.FREE +++ b/source/wham/src-M/DIMENSIONS.FREE @@ -3,7 +3,7 @@ integer MaxR,MaxT_h integer MaxSlice parameter (Max_Parm=1) - parameter (MaxQ=1,MaxQ1=MaxQ+2) + parameter (MaxQ=1,MaxQ1=MaxQ+3) parameter(MaxR=1,MaxT_h=32) parameter(MaxSlice=40) integer MaxN diff --git a/source/wham/src-M/enecalc1.F b/source/wham/src-M/enecalc1.F index c4f1f45..697c31d 100644 --- a/source/wham/src-M/enecalc1.F +++ b/source/wham/src-M/enecalc1.F @@ -91,6 +91,11 @@ q(nQ+1,iii+1)=rmsnat(iii+1) endif q(nQ+2,iii+1)=gyrate(iii+1) + if (q(nQ+1,iii+1).le.frac_cutoff) then + q(nQ+3,iii+1)=1.0d0 + else + q(nQ+3,iii+1)=0.0d0 + endif c fT=T0*beta_h(ib,ipar)*1.987D-3 c ft=2.0d0/(1.0d0+1.0d0/(T0*beta_h(ib,ipar)*1.987D-3)) if (rescale_mode.eq.1) then diff --git a/source/wham/src-M/readrtns.F b/source/wham/src-M/readrtns.F index 5dd092e..6762578 100644 --- a/source/wham/src-M/readrtns.F +++ b/source/wham/src-M/readrtns.F @@ -51,6 +51,8 @@ call multreadi(controlcard,"ISAMPL",isampl,nparmset,1) write (iout,*) "MaxSlice",MaxSlice call readi(controlcard,"NSLICE",nslice,1) + call reada(controlcard,"FRAC_CUTOFF",frac_cutoff,10.0d0) + write (iout,*) "RMSD cutoff for fraction calculations",frac_cutoff call flush(iout) if (nslice.gt.MaxSlice) then write (iout,*) "Error: parameter out of range: NSLICE",nslice, diff --git a/source/wham/src-M/wham_calc1.F b/source/wham/src-M/wham_calc1.F index 6768dcf..ff78e8d 100644 --- a/source/wham/src-M/wham_calc1.F +++ b/source/wham/src-M/wham_calc1.F @@ -531,7 +531,7 @@ c write (iout,*) 1.0d0/(beta_h(ib,iparm)*1.987D-3),ft sumE_p(i,iparm)=0.0d0 sumEbis_p(i,iparm)=0.0d0 sumEsq_p(i,iparm)=0.0d0 - do j=1,nQ+2 + do j=1,nQ+3 sumQ_p(j,i,iparm)=0.0d0 sumQsq_p(j,i,iparm)=0.0d0 sumEQ_p(j,i,iparm)=0.0d0 @@ -546,7 +546,7 @@ c write (iout,*) 1.0d0/(beta_h(ib,iparm)*1.987D-3),ft sumE(i,iparm)=0.0d0 sumEbis(i,iparm)=0.0d0 sumEsq(i,iparm)=0.0d0 - do j=1,nQ+2 + do j=1,nQ+3 sumQ(j,i,iparm)=0.0d0 sumQsq(j,i,iparm)=0.0d0 sumEQ(j,i,iparm)=0.0d0 @@ -826,7 +826,7 @@ c write (iout,*) "ftbis",ftbis sumE_p(k,iparm)=sumE_p(k,iparm)+etot*weight sumEbis_p(k,iparm)=sumEbis_p(k,iparm)+ebis*weight sumEsq_p(k,iparm)=sumEsq_p(k,iparm)+etot**2*weight - do j=1,nQ+2 + do j=1,nQ+3 sumQ_p(j,k,iparm)=sumQ_p(j,k,iparm)+q(j,t)*weight sumQsq_p(j,k,iparm)=sumQsq_p(j,k,iparm)+q(j,t)**2*weight sumEQ_p(j,k,iparm)=sumEQ_p(j,k,iparm) @@ -837,7 +837,7 @@ c write (iout,*) "ftbis",ftbis sumE(k,iparm)=sumE(k,iparm)+etot*weight sumEbis(k,iparm)=sumEbis(k,iparm)+ebis*weight sumEsq(k,iparm)=sumEsq(k,iparm)+etot**2*weight - do j=1,nQ+2 + do j=1,nQ+3 sumQ(j,k,iparm)=sumQ(j,k,iparm)+q(j,t)*weight sumQsq(j,k,iparm)=sumQsq(j,k,iparm)+q(j,t)**2*weight sumEQ(j,k,iparm)=sumEQ(j,k,iparm) @@ -1087,7 +1087,7 @@ c write (iout,*) "ftbis",ftbis & sumW(i,iparm) sumEsq(i,iparm)=(sumEsq(i,iparm)/sumW(i,iparm) & -sumE(i,iparm)**2)/(1.987D-3*(startGridT+i*delta_T)**2) - do j=1,nQ+2 + do j=1,nQ+3 sumQ(j,i,iparm)=sumQ(j,i,iparm)/sumW(i,iparm) sumQsq(j,i,iparm)=sumQsq(j,i,iparm)/sumW(i,iparm) & -sumQ(j,i,iparm)**2 @@ -1098,15 +1098,15 @@ c write (iout,*) "ftbis",ftbis & (startGridT+i*delta_T))+potEmin write (iout,'(f7.1,2f15.5,$)') startGridT+i*delta_T, & sumW(i,iparm),sumE(i,iparm) - write (iout,'(f10.5,$)') (sumQ(j,i,iparm),j=1,nQ+2) + write (iout,'(f10.5,$)') (sumQ(j,i,iparm),j=1,nQ+3) write (iout,'(e15.5,$)') sumEsq(i,iparm)-sumEbis(i,iparm), - & (sumQsq(j,i,iparm),j=1,nQ+2),(sumEQ(j,i,iparm),j=1,nQ+2) + & (sumQsq(j,i,iparm),j=1,nQ+3),(sumEQ(j,i,iparm),j=1,nQ+3) write (iout,*) write (34,'(f7.1,2f15.5,$)') startGridT+i*delta_T, & sumW(i,iparm),sumE(i,iparm) - write (34,'(f10.5,$)') (sumQ(j,i,iparm),j=1,nQ+2) + write (34,'(f10.5,$)') (sumQ(j,i,iparm),j=1,nQ+3) write (34,'(e15.5,$)') sumEsq(i,iparm)-sumEbis(i,iparm), - & (sumQsq(j,i,iparm),j=1,nQ+2),(sumEQ(j,i,iparm),j=1,nQ+2) + & (sumQsq(j,i,iparm),j=1,nQ+3),(sumEQ(j,i,iparm),j=1,nQ+3) write (34,*) call flush(34) enddo