X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Funres%2Fsrc_Eshel%2Fdihed_cons.F;fp=source%2Funres%2Fsrc_Eshel%2Fdihed_cons.F;h=e45405f28565bf2c40d7a5b2f8f0b3e0c366f8a4;hb=7308760ff07636ef6b1ee28d8c3a67a23c14b34b;hp=0000000000000000000000000000000000000000;hpb=9a54ab407f6d0d9d564d52763b3e2136450b9ffc;p=unres.git

diff --git a/source/unres/src_Eshel/dihed_cons.F b/source/unres/src_Eshel/dihed_cons.F
new file mode 100644
index 0000000..e45405f
--- /dev/null
+++ b/source/unres/src_Eshel/dihed_cons.F
@@ -0,0 +1,185 @@
+      subroutine secstrp2dihc
+      implicit real*8 (a-h,o-z)
+      include 'DIMENSIONS'
+      include 'COMMON.GEO'
+      include 'COMMON.BOUNDS'
+      include 'COMMON.CHAIN'
+      include 'COMMON.TORCNSTR'
+	include 'COMMON.IOUNITS'
+      character*1 secstruc(maxres)
+      COMMON/SECONDARYS/secstruc
+      character*80 line
+      logical errflag
+      external ilen
+
+cdr      call getenv_loc('SECPREDFIL',secpred)
+      lenpre=ilen(prefix)
+      secpred=prefix(:lenpre)//'.spred'
+
+#if defined(WINIFL) || defined(WINPGI)
+      open(isecpred,file=secpred,status='old',readonly,shared)
+#elif (defined CRAY) || (defined AIX)
+      open(isecpred,file=secpred,status='old',action='read')
+#elif (defined G77)
+      open(isecpred,file=secpred,status='old')
+#else
+      open(isecpred,file=secpred,status='old',action='read')
+#endif
+C read secondary structure prediction from JPRED here!
+!      read(isecpred,'(A80)',err=100,end=100) line
+!      read(line,'(f10.3)',err=110) ftors
+       read(isecpred,'(f10.3)',err=110) ftors
+
+      write (iout,*) 'FTORS factor =',ftors
+! initialize secstruc to any
+       do i=1,nres
+        secstruc(i) ='-'
+      enddo
+      ndih_constr=0
+      ndih_nconstr=0
+
+      call read_secstr_pred(isecpred,iout,errflag)
+      if (errflag) then
+         write(iout,*)'There is a problem with the list of secondary-',
+     &     'structure prediction'
+         goto 100
+      endif
+C 8/13/98 Set limits to generating the dihedral angles
+      do i=1,nres
+        phibound(1,i)=-pi
+        phibound(2,i)=pi
+      enddo
+      
+      ii=0
+      do i=1,nres
+        if ( secstruc(i) .eq. 'H') then
+C Helix restraints for this residue               
+           ii=ii+1 
+           idih_constr(ii)=i
+           phi0(ii) = 45.0D0*deg2rad
+           drange(ii)= 5.0D0*deg2rad
+           phibound(1,i) = phi0(ii)-drange(ii)
+           phibound(2,i) = phi0(ii)+drange(ii)
+        else if (secstruc(i) .eq. 'E') then
+C strand restraints for this residue               
+           ii=ii+1 
+           idih_constr(ii)=i 
+           phi0(ii) = 180.0D0*deg2rad
+           drange(ii)= 5.0D0*deg2rad
+           phibound(1,i) = phi0(ii)-drange(ii)
+           phibound(2,i) = phi0(ii)+drange(ii)
+        else
+C no restraints for this residue               
+           ndih_nconstr=ndih_nconstr+1
+           idih_nconstr(ndih_nconstr)=i
+        endif        
+      enddo
+      ndih_constr=ii
+      return
+100   continue
+      write(iout,'(A30,A80)')'Error reading file SECPRED',secpred
+      return
+ 110  continue
+      write(iout,'(A20)')'Error reading FTORS'
+      return
+      end 
+
+      subroutine read_secstr_pred(jin,jout,errors)
+
+      implicit real*8 (a-h,o-z)
+      INCLUDE 'DIMENSIONS'
+      include 'COMMON.IOUNITS'
+      include 'COMMON.CHAIN'
+      character*1 secstruc(maxres)
+      COMMON/SECONDARYS/secstruc
+      EXTERNAL ILEN
+      character*80 line,line1,ucase
+      logical errflag,errors,blankline
+
+      errors=.false.
+      read (jin,'(a)') line
+      write (jout,'(2a)') '> ',line(1:78)
+      line1=ucase(line)
+C Remember that we number full residues starting from 2, then, iseq=1 and iseq=nres
+C correspond to the end-groups.  ADD to the secondary structure prediction "-" for the
+C end-groups in the input file "*.spred"
+
+      iseq=1
+      do while (index(line1,'$END').eq.0)
+* Override commented lines.
+         ipos=1
+         blankline=.false.
+         do while (.not.blankline)
+            line1=' '
+            call mykey(line,line1,ipos,blankline,errflag) 
+            if (errflag) write (jout,'(2a)')
+     & 'Error when reading sequence in line: ',line
+            errors=errors .or. errflag
+            if (.not. blankline .and. .not. errflag) then
+               ipos1=2
+               iend=ilen(line1)
+               if (iseq.le.maxres) then
+                  if (line1(1:1).eq.'-' ) then 
+                     secstruc(iseq)=line1(1:1)
+                  else if ( ( ucase(line1(1:1)).eq.'E' ) .or.
+     &                      ( ucase(line1(1:1)).eq.'H' ) ) then
+                     secstruc(iseq)=ucase(line1(1:1))
+                  else
+                     errors=.true.
+                     write (jout,1010) line1(1:1), iseq
+                     goto 80
+                  endif                  
+               else
+                  errors=.true.
+                  write (jout,1000) iseq,maxres
+                  goto 80
+               endif
+               do while (ipos1.le.iend)
+
+                  iseq=iseq+1
+                  il=1
+                  ipos1=ipos1+1
+                  if (iseq.le.maxres) then
+                     if (line1(ipos1-1:ipos1-1).eq.'-' ) then 
+                        secstruc(iseq)=line1(ipos1-1:ipos1-1)
+                     else if((ucase(line1(ipos1-1:ipos1-1)).eq.'E').or.
+     &                     (ucase(line1(ipos1-1:ipos1-1)).eq.'H') ) then
+                        secstruc(iseq)=ucase(line1(ipos1-1:ipos1-1))
+                     else
+                        errors=.true.
+                        write (jout,1010) line1(ipos1-1:ipos1-1), iseq
+                        goto 80
+                     endif                  
+                  else
+                     errors=.true.
+                     write (jout,1000) iseq,maxres
+                     goto 80
+                  endif
+               enddo
+               iseq=iseq+1
+            endif
+         enddo
+         read (jin,'(a)') line
+         write (jout,'(2a)') '> ',line(1:78)
+         line1=ucase(line)
+      enddo
+
+cd    write (jout,'(10a8)') (sequence(i),i=1,iseq-1)
+
+cd check whether the found length of the chain is correct.
+      length_of_chain=iseq-1
+      if (length_of_chain .ne. nres) then
+!        errors=.true.
+        write (jout,'(a,i4,a,i4,a)') 
+     & 'Error: the number of labels specified in $SEC_STRUC_PRED ('
+     & ,length_of_chain,') does not match with the number of residues ('
+     & ,nres,').' 
+      endif
+   80 continue
+
+ 1000 format('Error - the number of residues (',i4,
+     & ') has exceeded maximum (',i4,').')
+ 1010 format ('Error - unrecognized secondary structure label',a4,
+     & ' in position',i4)
+      return
+      end