Cluster/unres modified to read MD coordinates in CX format, also from multiple files.
authorAdam Liwo <adam@piasek4.chem.univ.gda.pl>
Tue, 11 Aug 2015 04:09:06 +0000 (06:09 +0200)
committerAdam Liwo <adam@piasek4.chem.univ.gda.pl>
Tue, 11 Aug 2015 04:09:06 +0000 (06:09 +0200)
14 files changed:
bin/cluster/unres_clust-ifort.exe
source/cluster/unres/src/COMMON.CONTROL
source/cluster/unres/src/COMMON.IOUNITS
source/cluster/unres/src/DIMENSIONS
source/cluster/unres/src/Makefile_ifort
source/cluster/unres/src/compinfo.c
source/cluster/unres/src/fitsq.f
source/cluster/unres/src/hc.f
source/cluster/unres/src/intcor.f
source/cluster/unres/src/main_clust.f [deleted file]
source/cluster/unres/src/readrtns.F
source/cluster/unres/src/sizesclu.dat
source/cluster/unres/src/timing.F
source/cluster/unres/src/wrtclust.f

index f7afe5b..82bb05e 100755 (executable)
Binary files a/bin/cluster/unres_clust-ifort.exe and b/bin/cluster/unres_clust-ifort.exe differ
index 1fa1444..a9bf4e7 100644 (file)
@@ -1,7 +1,7 @@
-      double precision betaT
-      integer iscode,indpdb,outpdb,outmol2,iopt
+      double precision betaT,ts,te
+      integer iscode,indpdb,outpdb,outmol2,iopt,is,ie,isampl
       logical refstr,pdbref,punch_dist,print_dist,caonly,lside,
-     & lprint_cart,lprint_int,from_cart,efree
-      common /cntrl/ betaT,iscode,indpdb,refstr,pdbref,outpdb,outmol2,
-     & punch_dist,print_dist,caonly,lside,lprint_cart,lprint_int,
-     & from_cart,efree,iopt
+     & lprint_cart,lprint_int,from_cart,from_cx,efree
+      common /cntrl/ betaT,ts,te,is,ie,isampl,iscode,indpdb,refstr,
+     & pdbref,outpdb,outmol2,punch_dist,print_dist,caonly,lside,
+     & lprint_cart,lprint_int,from_cart,from_cx,efree,iopt
index b80aa37..334014f 100644 (file)
@@ -2,13 +2,14 @@ C-----------------------------------------------------------------------
 C I/O units used by the program
 C-----------------------------------------------------------------------
       integer inp,iout,igeom,intin,ipdb,imol2,ipdbin,jplot,jstatin,
-     & jstatout
+     & jstatout,nfiles
       common /iounits/ inp,iout,igeom,intin,ipdb,imol2,ipdbin,jplot,
-     & jrms,jstatin,jstatout
+     & jrms,jstatin,jstatout,nfiles
       character*256 outname,intname,intinname,pdbname,mol2name,prefix,
-     &  prefout,prefintin,statname,rmsname,statinname,statoutname
+     &  prefout,prefintin,statname,rmsname,statinname,statoutname,
+     &  cxfiles(maxfiles)
       common /fnames/ outname,intname,pdbname,mol2name,intinname,prefix,
-     &  prefout,prefintin,statinname,statoutname
+     &  prefout,prefintin,statinname,statoutname,cxfiles
 C-----------------------------------------------------------------------
 C INP    - main input file
 C IOUT   - list file
index 44a47a7..8b077ff 100644 (file)
@@ -71,3 +71,5 @@ C Maximum number of chains
       parameter (mxch=1)
 C Maximum number of generated conformations
       parameter (mxio=1000)
+C Maximum number of coordinate files
+      parameter (maxfiles=10)
index 19e93e8..595c461 100644 (file)
@@ -1,16 +1,18 @@
 BIN=/home/adam/UNRES/bin
 INSTALL_DIR = /usr/local/mpich-1.2.7p1-intel
-FC= ifort
-OPT = -O3 -ip -w 
-FFLAGS =  ${OPT} -I. 
-CPPFLAGS = -DLINUX
+#FC= ifort
+#OPT = -O3 -ip -w 
+#OPT = -g -CB
+#FFLAGS =  ${OPT} -I. 
+#CPPFLAGS = -DLINUX -DNEWUNRES
 
 # Next for debugging
-#FC = f77
+FC = gfortran
 #OPT = -g -fbounds-check
-#CPPFLAGS = -DLINUX -DG77
+OPT = -O3
+CPPFLAGS = -DLINUX -DG77 -DNEWUNRES
 
-LIBS = 
+LIBS = xdrf/libxdrf.a
 
 .SUFFIXES: .F
 .F.o:
@@ -18,16 +20,19 @@ LIBS =
 .f.o:
        $(FC) -c $(CPPFLAGS) $(FFLAGS) $<
 
-objects = main_clust.o arcos.o cartprint.o chainbuild.o convert.o initialize.o \
+objects = main_clust.o arcos.o cartprint.o chainbuild.o convert.o initialize.o cxread.o \
        matmult.o readrtns.o pinorm.o rescode.o intcor.o timing.o misc.o \
        geomout.o readpdb.o fitsq.o hc.o  \
        track.o wrtclust.o srtclust.o noyes.o contact.o
 
-unres_clust: $(objects)
+unres_clust: $(objects) xdrf/libxdrf.a
        cc -o compinfo compinfo.c
        ./compinfo | true
        ${FC} -c cinfo.f
-       $(FC) ${OPT} ${objects} cinfo.o ${LIBS} -o ../../../../bin/cluster/unres_clust-ifort.exe
+       $(FC) ${OPT} ${objects} cinfo.o ${LIBS} -o ../../../../bin/cluster/unres_clust-ifort-test.exe
+
+xdrf/libxdrf.a:
+       cd xdrf && make
 
 clean:
        /bin/rm *.o
index e28f686..af53dc0 100644 (file)
@@ -21,6 +21,7 @@ if (fgets(buf,498,in) != NULL)
 iv3++;
 fprintf(out,"C %d %d %d\n",iv1,iv2,iv3);
 fprintf(out,"      subroutine cinfo\n");
+fprintf(out,"      include 'DIMENSIONS'\n");
 fprintf(out,"      include 'COMMON.IOUNITS'\n");
 fprintf(out,"      write(iout,*)'++++ Compile info ++++'\n");
 fprintf(out,"      write(iout,*)'Version %d.%-d build %d'\n",iv1,iv2,iv3);
index 17d92ee..fbcf223 100644 (file)
@@ -1,5 +1,5 @@
       subroutine fitsq(rms,x,y,nn,t,b,non_conv)
-      implicit real*8 (a-h,o-z)
+      include 'DIMENSIONS'
       include 'COMMON.IOUNITS'
 c  x and y are the vectors of coordinates (dimensioned (3,n)) of the two
 c  structures to be superimposed.  nn is 3*n, where n is the number of  
index 2e03f24..da10b8f 100644 (file)
@@ -282,6 +282,7 @@ C statement was not necessary.                                  C
 C---------------------------------------------------------------C
       SUBROUTINE HCASS(N,IA,IB,CRIT,LEV,ICLASS,HVALS,IORDER,
      X        CRITVAL,HEIGHT)
+      include 'DIMENSIONS'
       include 'sizesclu.dat'
       include 'COMMON.IOUNITS'
       integer ICLASS(maxconf,maxconf-1)
@@ -410,6 +411,7 @@ C  F. Murtagh, ESA/ESO/STECF, Garching, Feb. 1986.C
 C                                                 C 
 C-------------------------------------------------C
       SUBROUTINE HCDEN(LEV,IORDER,HEIGHT,CRITVAL)
+      include 'DIMENSIONS'
       include 'COMMON.IOUNITS'
       CHARACTER*80 LINE
       INTEGER IORDER(LEV),HEIGHT(LEV)
index 7d0b849..045cce1 100644 (file)
@@ -8,7 +8,7 @@ c
       include 'DIMENSIONS'
       include 'COMMON.GEO'
       include 'COMMON.CHAIN'
-      if (i1.lt.1 .or. i2.lt.1 .or. i3.lt.1 .or. i4.lt.1) then
+      if (i1.lt.1 .or. i2.lt.1 .or. i3.lt.1) then
         alpha=0.0
         return
       endif
diff --git a/source/cluster/unres/src/main_clust.f b/source/cluster/unres/src/main_clust.f
deleted file mode 100644 (file)
index e3dd395..0000000
+++ /dev/null
@@ -1,498 +0,0 @@
-C
-C Program to cluster united-residue MCM results.
-C
-      include 'DIMENSIONS'
-      include 'sizesclu.dat'
-      include 'COMMON.TIME1'
-      include 'COMMON.INTERACT'
-      include 'COMMON.NAMES'
-      include 'COMMON.GEO'
-      include 'COMMON.HEADER'
-      include 'COMMON.CONTROL'
-      include 'COMMON.CONTACTS'
-      include 'COMMON.CHAIN'
-      include 'COMMON.VAR'
-      include 'COMMON.CLUSTER'
-      include 'COMMON.IOUNITS'
-      logical printang(max_cut)
-      integer printpdb(max_cut)
-      integer printmol2(max_cut)
-      character*240 lineh
-      REAL CRIT(maxconf),MEMBR(maxconf)
-      REAL CRITVAL(maxconf-1)
-      REAL DISS(maxdist)
-      INTEGER IA(maxconf),IB(maxconf)
-      INTEGER ICLASS(maxconf,maxconf-1),HVALS(maxconf-1)
-      INTEGER IORDER(maxconf-1),HEIGHT(maxconf-1)
-      DIMENSION NN(maxconf),DISNN(maxconf)
-      LOGICAL FLAG(maxconf)
-      
-      double precision varia(maxvar)
-      double precision hrtime,mintime,sectime
-      logical eof
-      double precision eee(n_ene,maxconf)
-
-      call initialize
-      call openunits
-      call read_control
-      call molread
-      do i=1,nres
-        phi(i)=0.0D0
-        theta(i)=0.0D0
-        alph(i)=0.0D0
-        omeg(i)=0.0D0
-      enddo
-
-      print *,'MAIN: nnt=',nnt,' nct=',nct
-
-      DO I=1,NCUT
-        PRINTANG(I)=.FALSE.
-        PRINTPDB(I)=0
-        printmol2(i)=0
-        IF (RCUTOFF(I).LT.0.0) THEN
-          RCUTOFF(I)=ABS(RCUTOFF(I))
-          PRINTANG(I)=.TRUE.
-          PRINTPDB(I)=outpdb
-          printmol2(i)=outmol2
-        ENDIF
-      ENDDO
-      PRINT *,'Number of cutoffs:',NCUT
-      PRINT *,'Cutoff values:'
-      DO ICUT=1,NCUT
-        PRINT '(8HRCUTOFF(,I2,2H)=,F8.1,2i2)',ICUT,RCUTOFF(ICUT),
-     &    printpdb(icut),printmol2(icut)
-      ENDDO
-      DO I=1,NRES-3  
-        MULT(I)=1
-      ENDDO
-      ICON=1
-  123 continue
-      if (from_cart) then
-        if (efree) then
-        read (intin,*,end=13,err=11) energy(icon),totfree(icon),
-     &    rmstab(icon),
-     &    nss_all(icon),(ihpb_all(ii,icon),jhpb_all(i,icon),
-     &    i=1,nss_all(icon)),iscore(icon)
-        else
-        read (intin,*,end=13,err=11) energy(icon),rmstab(icon),
-     &    nss_all(icon),(ihpb_all(ii,icon),jhpb_all(i,icon),
-     &    i=1,nss_all(icon)),iscore(icon)
-        endif
-        read (intin,'(8f10.5)',end=13,err=10) 
-     &    ((allcart(j,i,icon),j=1,3),i=1,nres),
-     &    ((allcart(j,i+nres,icon),j=1,3),i=nnt,nct)
-        print *,icon,energy(icon),nss_all(icon),rmstab(icon)
-      else 
-        read(intin,'(a80)',end=13,err=12) lineh
-        read(lineh(:5),*,err=8) ic
-        if (efree) then
-        read(lineh(6:),*,err=8) energy(icon)
-        else
-        read(lineh(6:),*,err=8) energy(icon)
-        endif
-        goto 9
-    8   ic=1
-        print *,'error, assuming e=1d10',lineh
-        energy(icon)=1d10
-        nss=0
-    9   continue
-cold        read(lineh(18:),*,end=13,err=11) nss_all(icon)
-        ii = index(lineh(15:)," ")+15
-        read(lineh(ii:),*,end=13,err=11) nss_all(icon)
-        IF (NSS_all(icon).LT.9) THEN
-          read (lineh(20:),*,end=102)
-     &    (IHPB_all(I,icon),JHPB_all(I,icon),I=1,NSS_all(icon)),
-     &    iscore(icon)
-        ELSE
-          read (lineh(20:),*,end=102) 
-     &           (IHPB_all(I,icon),JHPB_all(I,icon),I=1,8)
-          read (intin,*) (IHPB_all(I,icon),JHPB_all(I,icon),
-     &      I=9,NSS_all(icon)),iscore(icon)
-        ENDIF
-
-  102   continue  
-
-        PRINT *,'IC:',IC,' ENERGY:',ENERGY(ICON)
-        totfree(icon)=energy(icon)
-        call read_angles(intin,*13)
-        call chainbuild
-        do ii=1,2*nres
-          do jj=1,3
-            allcart(jj,ii,icon)=c(jj,ii)
-          enddo 
-        enddo
-        do i=1,nres
-          phiall(i,icon)=phi(i)
-          thetall(i,icon)=theta(i)
-          alphall(i,icon)=alph(i)
-          omall(i,icon)=omeg(i)
-        enddo
-      endif
-      ICON=ICON+1
-      GOTO 123
-C      
-C CALCULATE DISTANCES
-C
-   10 print *,'something wrong with angles'
-      goto 13
-   11 print *,'something wrong with NSS',nss
-      goto 13
-   12 print *,'something wrong with header'
-
-   13 NCON=ICON-1
-      if (lgrp) then
-        i=0
-        do while (i.lt.ncon)
-c          read (jstatin,*,err=1111,end=1111) idum,(eee(j,i),j=1,n_ene)
-          read (jstatin,'(a)') lineh
-          if (index(lineh,'#').gt.0) goto 1123
-          i=i+1
-          read (lineh,*) idum,(eee(j,i),j=1,n_ene)
-          print *,i,idum,(eee(j,i),j=1,n_ene)
-          energy(i)=eee(15,i) 
- 1123     continue
-        enddo
-      endif 
-      goto 1112
- 1111 print *,'Error in the statin file; statout will not be produced.'
-      write (iout,*)
-     &  'Error in the statin file; statout will not be produced.'
-      lgrp=.false.
- 1112 continue
-      print *,'ncon',ncon
-      DO I=1,NCON
-        ICC(I)=I
-      ENDDO
-      WRITE (iout,'(A80)') TITEL
-      t1=tcpu()
-      DO I=1,NCON-1
-        if (mod(i,100).eq.0) print *,'Calculating RMS i=',i
-        if (from_cart) then
-          do ii=1,2*nres
-            do jj=1,3
-              c(jj,ii)=allcart(jj,ii,i)
-            enddo 
-          enddo
-          call int_from_cart(.true.,.false.)
-          do ii=1,nres
-            phiall(ii,i)=phi(ii)
-            thetall(ii,i)=theta(ii)
-            alphall(ii,i)=alph(ii)
-            omall(ii,i)=omeg(ii)
-          enddo
-        else
-          do ii=1,nres
-            phi(ii)=phiall(ii,i)
-            theta(ii)=thetall(ii,i)
-            alph(ii)=alphall(ii,i)
-            omeg(ii)=omall(ii,i)
-          enddo
-          call chainbuild
-        endif
-        do ii=1,nres
-          do jj=1,3
-            cref(jj,ii)=c(jj,ii)
-          enddo
-        enddo
-        DO J=I+1,NCON
-          IND=IOFFSET(NCON,I,J)
-          ATTALUMS(IND)=DIFCONF(I,J)
-c          write (iout,'(2i4,f10.5)') i,j,ATTALUMS(IND)
-          DISS(IND)=ATTALUMS(IND)
-        ENDDO
-      ENDDO
-      t2=tcpu()
-      WRITE (iout,'(/a,1pe14.5,a/)') 
-     & 'Time for distance calculation:',T2-T1,' sec.'
-      t1=tcpu()
-      PRINT '(a)','End of distance computation'
-
-      if (punch_dist) then
-        do i=1,ncon-1
-          do j=i+1,ncon
-            IND=IOFFSET(NCON,I,J)
-            write (jrms,'(2i5,2f10.5)') i,j,attalums(IND),
-     &        energy(j)-energy(i)
-          enddo
-        enddo
-      endif
-C
-C Print out the RMS deviation matrix.
-C
-      if (print_dist) CALL DISTOUT(NCON)
-C
-C  call hierarchical clustering HC from F. Murtagh
-C
-      N=NCON
-      LEN = (N*(N-1))/2
-      write(iout,*) "-------------------------------------------"
-      write(iout,*) "HIERARCHICAL CLUSTERING using"
-      if (iopt.eq.1) then
-        write(iout,*) "WARD'S MINIMUM VARIANCE METHOD"
-      elseif (iopt.eq.2) then
-        write(iout,*) "SINGLE LINK METHOD"
-      elseif (iopt.eq.3) then
-        write(iout,*) "COMPLETE LINK METHOD"
-      elseif (iopt.eq.4) then
-        write(iout,*) "AVERAGE LINK (OR GROUP AVERAGE) METHOD"
-      elseif (iopt.eq.5) then
-        write(iout,*) "MCQUITTY'S METHOD"
-      elseif (iopt.eq.6) then
-        write(iout,*) "MEDIAN (GOWER'S) METHOD"
-      elseif (iopt.eq.7) then
-        write(iout,*) "CENTROID METHOD"
-      else
-        write(iout,*) "IOPT=",iopt," IS INVALID, use 1-7"
-        write(*,*) "IOPT=",iopt," IS INVALID, use 1-7"
-        stop
-      endif
-      write(iout,*)
-      write(iout,*) "hc.f by F. Murtagh, ESA/ESO/STECF, Garching"
-      write(iout,*) "February 1986"
-      write(iout,*) "References:"
-      write(iout,*) "1. Multidimensional clustering algorithms"
-      write(iout,*) "   Fionn Murtagh"
-      write(iout,*) "   Vienna : Physica-Verlag, 1985."
-      write(iout,*) "2. Multivariate data analysis"
-      write(iout,*) "   Fionn Murtagh and Andre Heck"
-      write(iout,*) "   Kluwer Academic Publishers, 1987"
-      write(iout,*) "-------------------------------------------"
-      write(iout,*)
-
-      CALL HC(N,M,LEN,IOPT,IA,IB,CRIT,MEMBR,NN,DISNN,FLAG,DISS)
-      LEV = N-1
-      CALL HCASS(N,IA,IB,CRIT,LEV,ICLASS,HVALS,IORDER,CRITVAL,HEIGHT)
-      CALL HCDEN(LEV,IORDER,HEIGHT,CRITVAL)
-
-      icut=1
-      i=1
-      NGR=i+1
-      do j=1,n
-        licz(iclass(j,i))=licz(iclass(j,i))+1
-        nconf(iclass(j,i),licz(iclass(j,i)))=j
-      enddo        
-      do i=1,lev-1
-
-         idum=lev-i
-         DO L=1,LEV
-            IF (HEIGHT(L).EQ.IDUM) GOTO 190
-         ENDDO
- 190     IDUM=L
-cd         print *,i+1,CRITVAL(IDUM)
-         IF (CRITVAL(IDUM).LT. RCUTOFF(ICUT)) THEN
-          WRITE (iout,'(/a,f10.5)') 'AT CUTOFF:',rcutoff(icut)
-          write (iout,'(a,f8.2)') 'Maximum distance found:',
-     &              CRITVAL(IDUM)
-          CALL SRTCLUST(ICUT,ncon)
-          CALL TRACK(ICUT)
-          CALL WRTCLUST(ncon,ICUT,PRINTANG,PRINTPDB,PRINTMOL2)
-          icut=icut+1
-          if (icut.gt.ncut) goto 191
-         ENDIF
-         NGR=i+1
-         do l=1,maxgr
-          licz(l)=0
-         enddo
-         do j=1,n
-          licz(iclass(j,i))=licz(iclass(j,i))+1
-          nconf(iclass(j,i),licz(iclass(j,i)))=j
-cd          print *,j,iclass(j,i),
-cd     &     licz(iclass(j,i)),nconf(iclass(j,i),licz(iclass(j,i)))
-         enddo
-      enddo
- 191  continue
-C
-      if (plot_tree) then
-        CALL WRITRACK
-        CALL PLOTREE
-      endif
-C
-      t2=tcpu()
-      WRITE (iout,'(/a,1pe14.5,a/)') 
-     & 'Total time for clustering:',T2-T1,' sec.'
-      if (lgrp) then
-        write (jstatout,'(a5,18a12,10f4.1)')"#    ",
-     &   "EVDW SC-SC","EVDW2 SC-p","EES p-p",
-     &   "EBE bending","ESC SCloc","ETORS ",
-     &   "ECORR4 ","ECORR5 ","ECORR6 ",
-     &   "EELLO ","ETURN3 ","ETURN4 ","ETURN6 "
-     &   ,"ETOT total","RMSD","nat.contact","nnt.contact",
-     &   "cont.order",(rcutoff(i),i=1,ncut)
-        do i=1,ncon
-          write(jstatout,'(i5,18f12.4,10i6)') i,(eee(j,i),j=1,n_ene),
-     &     (iass_tot(i,icut),icut=1,ncut)
-        enddo
-      endif
-C
-      stop '********** Program terminated normally.'
-      end
-c---------------------------------------------------------------------------
-      double precision function difconf(icon,jcon)
-      include 'DIMENSIONS'
-      include 'sizesclu.dat'
-      include 'COMMON.CONTROL'
-      include 'COMMON.CLUSTER'
-      include 'COMMON.CHAIN' 
-      include 'COMMON.INTERACT'
-      include 'COMMON.VAR'
-      logical non_conv
-      double precision przes(3),obrot(3,3)
-      double precision xx(3,maxres2),yy(3,maxres2)
-      do i=1,nres
-        phi(i)=phiall(i,jcon)
-        theta(i)=thetall(i,jcon)
-        alph(i)=alphall(i,jcon)
-        omeg(i)=omall(i,jcon)
-      enddo
-      call chainbuild
-c     do i=1,nres
-c       print '(i4,2(3f10.5,5x))',i,(cref(j,i),j=1,3),(c(j,i),j=1,3)
-c     enddo
-      if (lside) then
-        ii=0
-        do i=nnt,nct
-          ii=ii+1
-          do j=1,3
-            xx(j,ii)=allcart(j,i,jcon)
-            yy(j,ii)=cref(j,i)
-          enddo
-        enddo
-        do i=nnt,nct
-c          if (itype(i).ne.10) then
-            ii=ii+1
-            do j=1,3 
-              xx(j,ii)=allcart(j,i+nres,jcon)
-              yy(j,ii)=cref(j,i+nres)
-            enddo
-c          endif
-        enddo
-        call fitsq(rms,xx(1,1),yy(1,1),ii,przes,obrot,non_conv)
-      else
-        do i=nnt,nct
-          do j=1,3
-            c(j,i)=allcart(j,i,jcon)
-          enddo
-c          write (2,'(i5,3f10.5,5x,3f10.5)') i,
-c     &     (c(j,i),j=1,3),(cref(j,i),j=1,3)
-        enddo
-        call fitsq(rms,c(1,nnt),cref(1,nnt),nct-nnt+1,przes,obrot
-     &       ,non_conv)
-      endif
-      if (rms.lt.0.0) then
-        print *,error,'rms^2 = ',rms,icon,jcon
-        stop
-      endif
-      if (non_conv) print *,non_conv,icon,jcon
-      difconf=dsqrt(rms)
-      RETURN
-      END
-C------------------------------------------------------------------------------
-      double precision function rmsnat(jcon)
-      include 'DIMENSIONS'
-      include 'sizesclu.dat'
-      include 'COMMON.CONTROL'
-      include 'COMMON.CLUSTER'
-      include 'COMMON.CHAIN' 
-      include 'COMMON.INTERACT'
-      include 'COMMON.VAR'
-      logical non_conv
-      double precision przes(3),obrot(3,3)
-      double precision xx(3,maxres2),yy(3,maxres2)
-      if (lside) then
-        ii=0
-        do i=nnt,nct
-          ii=ii+1
-          do j=1,3
-            xx(j,ii)=allcart(j,i,jcon)
-            yy(j,ii)=cref_pdb(j,i)
-          enddo
-        enddo
-        do i=nnt,nct
-c          if (itype(i).ne.10) then
-            ii=ii+1
-            do j=1,3 
-              xx(j,ii)=allcart(j,i+nres,jcon)
-              yy(j,ii)=cref_pdb(j,i+nres)
-            enddo
-c          endif
-        enddo
-        call fitsq(rms,xx(1,1),yy(1,1),ii,przes,obrot,non_conv)
-      else
-        do i=nnt,nct
-          do j=1,3
-            c(j,i)=allcart(j,i,jcon)
-          enddo
-        enddo
-        call fitsq(rms,c(1,nnt),cref_pdb(1,nnt),nct-nnt+1,przes,obrot
-     &       ,non_conv)
-      endif
-      if (rms.lt.0.0) then
-        print *,error,'rms^2 = ',rms,icon,jcon
-        stop
-      endif
-      if (non_conv) print *,non_conv,icon,jcon
-      rmsnat=dsqrt(rms)
-      RETURN
-      END
-C------------------------------------------------------------------------------
-      subroutine distout(ncon)
-      include 'DIMENSIONS'
-      include 'sizesclu.dat'
-      parameter (ncol=10)
-      include 'COMMON.IOUNITS'
-      include 'COMMON.CLUSTER'
-      dimension b(ncol)
-      write (iout,'(a)') 'The distance matrix'
-      do 1 i=1,ncon,ncol
-      nlim=min0(i+ncol-1,ncon)
-      write (iout,1000) (k,k=i,nlim)
-      write (iout,'(8h--------,10a)') ('-------',k=i,nlim)
- 1000 format (/8x,10(i4,3x))
- 1020 format (/1x,80(1h-)/)
-      do 2 j=i,ncon
-      jlim=min0(j,nlim)
-      if (jlim.eq.j) then
-        b(jlim-i+1)=0.0d0
-        jlim1=jlim-1
-      else
-        jlim1=jlim
-      endif
-      do 3 k=i,jlim1
-       if (j.lt.k) then 
-          IND=IOFFSET(NCON,j,k)
-       else
-          IND=IOFFSET(NCON,k,j)
-       endif
-    3  b(k-i+1)=attalums(IND)
-      write (iout,1010) j,(b(k),k=1,jlim-i+1)
-    2 continue
-    1 continue
- 1010 format (i5,3x,10(f6.2,1x))
-      return
-      end
-C-----------------------------------------------------------------------
-      block data
-      include 'DIMENSIONS'
-      include 'COMMON.LOCAL'
-      data dsc /  1.237,  ! CYS (type  1)
-     &            2.142,  ! MET (type  2)
-     &            2.299,  ! PHE (type  3)
-     &            1.776,  ! ILE (type  4)
-     &            1.939,  ! LEU (type  5)
-     &            1.410,  ! VAL (type  6)
-     &            2.605,  ! TRP (type  7)
-     &            2.484,  ! TYR (type  8)
-     &            0.743,  ! ALA (type  9)
-     &            0.000,  ! GLY (type 10)
-     &            1.393,  ! THR (type 11)
-     &            1.150,  ! SER (type 12)
-     &            2.240,  ! GLN (type 13)
-     &            1.684,  ! ASN (type 14)
-     &            2.254,  ! GLU (type 15)
-     &            1.709,  ! ASP (type 16)
-     &            2.113,  ! HIS (type 17)
-     &            3.020,  ! ARG (type 18)
-     &            2.541,  ! LYS (type 19)
-     &            1.345,  ! PRO (type 20)
-     &            0.000  /! D   (type 21)
-      end
index 9d55005..86d6dae 100644 (file)
@@ -48,10 +48,22 @@ C Set up the time limit (caution! The time must be input in minutes!)
       lgrp=(index(controlcard,'LGRP').gt.0)
       caonly=(index(controlcard,'CA_ONLY').gt.0)
       print_dist=(index(controlcard,'PRINT_DIST').gt.0)
+      print *,"print_dist",print_dist
       call multreada(controlcard,'CUTOFF',rcutoff,ncut,-1.0d0)
       call readi(controlcard,'IOPT',iopt,2) 
+      print *,"iopt",iopt
       lside = index(controlcard,"SIDE").gt.0
       efree = index(controlcard,"EFREE").gt.0
+      print *,controlcard
+      call readi(controlcard,'ISTART',is,1)
+      print *,"is",is
+      call readi(controlcard,'IEND',ie,10000000)
+      print *,"ie",ie
+      call readi(controlcard,'ISAMPL',isampl,1)
+      print *,"isampl",isampl
+      call reada(controlcard,'TS',ts,1.0d10)
+      call reada(controlcard,'TE',te,1.0d10)
+      print *,"is",is," ie",ie
       if (min_var) iopt=1
       return
       end
@@ -353,6 +365,37 @@ c----------------------------------------------------------------------------
       return
       end
 c----------------------------------------------------------------------------
+      subroutine split_string(rekord,tablica,dim,nsub)
+      implicit none
+      integer dim,nsub,i,ii,ll,kk
+      character*(*) tablica(dim)
+      character*(*) rekord
+      integer ilen
+      external ilen
+      do i=1,dim
+        tablica(i)=" "
+      enddo
+      ii=1
+      ll = ilen(rekord)
+      nsub=0
+      do i=1,dim
+C Find the start of term name
+        kk = 0
+        do while (ii.le.ll .and. rekord(ii:ii).eq." ")
+          ii = ii+1
+        enddo
+C Parse the name into TABLICA(i) until blank found
+        do while (ii.le.ll .and. rekord(ii:ii).ne." ")
+          kk = kk+1
+          tablica(i)(kk:kk)=rekord(ii:ii)
+          ii = ii+1
+        enddo
+        if (kk.gt.0) nsub=nsub+1
+        if (ii.gt.ll) return
+      enddo
+      return
+      end
+c----------------------------------------------------------------------------
       subroutine card_concat(card)
       include 'DIMENSIONS'
       include 'COMMON.IOUNITS'
@@ -385,26 +428,43 @@ c----------------------------------------------------------------------------
       call getenv('PDB',cfrom_pdb)
       call getenv('PRINTCOOR',cprint)
       from_cart = index(ucase(cfrom_pdb),'CART').gt.0
-      lprint_cart = index(ucase(cprint),'PRINT_CART').gt.0
-      lprint_int = index(ucase(cprint),'PRINT_INT').gt.0
-      if (from_cart  .and. .not.lprint_int) then
-        lprint_cart=.true.
-        lprint_int=.false.
-      endif
-      if (.not.lprint_cart .and. .not.lprint_int) lprint_int=.true.
+      from_cx = index(ucase(cfrom_pdb),'CX').gt.0
+      lprint_cart = index(ucase(cprint),'CART').gt.0
+      lprint_int = index(ucase(cprint),'INT').gt.0
+c      if (.not.lprint_cart .and. .not.lprint_int) lprint_int=.true.
       lenpre=ilen(prefix)
       lenout=ilen(prefout)
       lenint=ilen(prefintin)
 C Get the names and open the input files
       open (inp,file=prefix(:ilen(prefix))//'.inp',status='old')
       outname=prefout(:lenout)//'_clust.out'
-      if (lprint_cart) then
-        intname=prefintin(:lenint)//'_clust'//'.x'
+      if (from_cart) then
         intinname=prefintin(:lenint)//'.x'
+      else if (from_cx) then
+        write (iout,*) "cx files: ",prefintin(:ilen(prefintin))
+        call split_string(prefintin,cxfiles(1),maxfiles,nfiles)
+        write (iout,*) "nfiles",nfiles
+        write (iout,*) "Split cxfiles"
+        do i=1,nfiles
+          cxfiles(i)=cxfiles(i)(:ilen(cxfiles(i)))//'.cx'
+          write (iout,*) cxfiles(i)(:ilen(cxfiles(i)))
+        enddo
       else
-        intname=prefintin(:lenint)//'_clust'//'.int'
         intinname=prefintin(:lenint)//'.int'
       endif
+      if (lprint_cart) then
+        if (.not. from_cx) then
+          intname=prefintin(:lenint)//'_clust'//'.x'
+        else
+          intname=cxfiles(1)(:ilen(cxfiles(1)))//'_clust'//'.x'
+        endif
+      else if (lprint_int) then
+        if (.not. from_cx) then
+          intname=prefintin(:lenint)//'_clust'//'.int'
+        else
+          intname=cxfiles(1)(:ilen(cxfiles(1)))//'_clust'//'.int'
+        endif
+      endif
       rmsname=prefintin(:lenint)//'.rms'
       statinname=prefintin(:lenint)//'.stat'
       print *,statinname
@@ -412,7 +472,7 @@ C Get the names and open the input files
       open (jplot,file=prefout(:ilen(prefout))//'.tex',
      &       status='unknown')
       print *,'unit',jplot,' opened'
-      open (intin,file=intinname,status='old')
+      if (.not. from_cx) open (intin,file=intinname,status='old')
       print *,'unit',intin,' opened'
       open (jrms,file=rmsname,status='unknown')
       open (jstatin,file=statinname,status='unknown')
index c58a5f2..cb6c733 100644 (file)
@@ -4,7 +4,7 @@
 *
 * Max. number of conformations in the data set.
 *
-      PARAMETER (MAXCONF=2500)
+      PARAMETER (MAXCONF=5000)
 *
 * Max. number of "distances" between conformations.
 *
index ecab9b4..9bca353 100644 (file)
@@ -170,6 +170,7 @@ C Next definitions for Linux
 C=========================================================================
 C
       subroutine dajczas(rntime,hrtime,mintime,sectime)
+      include 'DIMENSIONS'
       include 'COMMON.IOUNITS'
       real*8 rntime,hrtime,mintime,sectime 
       hrtime=rntime/3600.0D0 
index ca43876..3f057b2 100644 (file)
@@ -25,20 +25,20 @@ c      ICANT(I,J)=((NCON+NCON-J)*(J-1))/2+I-J
 C
 C  PRINT OUT THE RESULTS OF CLUSTER ANALYSIS
 C
-      ii1= index(intinname,'/')
+      ii1= index(intname,'/')
       ii2=ii1
       ii1=ii1+1
       do while (ii2.gt.0) 
         ii1=ii1+ii2
-        ii2=index(intinname(ii1:),'/')
+        ii2=index(intname(ii1:),'/')
       enddo 
-      ii = ii1+index(intinname(ii1:),'.')-1
+      ii = ii1+index(intname(ii1:),'.')-1
       if (ii.eq.0) then
-        ii=ilen(intinname)
+        ii=ilen(intname)
       else
         ii=ii-1
       endif
-      prefixp=intinname(ii1:ii)
+      prefixp=intname(ii1:ii)
 cd    print *,icut,printang(icut),printpdb(icut),printmol2(icut)
 cd    print *,'ecut=',ecut
       WRITE (iout,100) NGR
@@ -94,9 +94,17 @@ cd          print '(3i4,f12.4)',ind,ii,jj,curr_dist
       write (iout,'(/A,F8.1,A,F8.1)')
      & 'Max. distance in the family:',amax_dim,
      & '; average distance in the family:',ave_dim 
-      if (refstr .or. pdbref) write (iout,'(a,i5,f8.3)') 
+      if (refstr .or. pdbref) then 
+        write (iout,'(a,i5,2f8.3)')
      & "RMSD of the lowest-energy conformation #",nconf(igr,1),
-     &  rmsnat(nconf(igr,1))
+     &  rmsnat(nconf(igr,1)),rmstab(nconf(igr,1))
+        rmsave=0.0d0
+        do i=1,licz(igr)
+          rmsave=rmsave+rmsnat(nconf(igr,i))
+        enddo
+        rmsave=rmsave/licz(igr)
+        write (iout,'(a,f8.3)') "Average RMSD in the family",rmsave
+      endif
    19 CONTINUE
       WRITE (iout,400)
       WRITE (iout,500) (I,IASS(I),I=1,NCON)
@@ -237,6 +245,16 @@ c Write conformations of the family i to PDB files
           else
 c Produce only a single PDB file for the leading member of the family
             write (iout,*) 'Writing pdb file: icon=',icon
+            if (from_cart .or. from_cx) then
+            
+            do ii=1,2*nres
+              do j=1,3
+              c(j,ii)=allcart(j,ii,icon)
+              enddo
+            enddo
+
+            else
+
             do ii=1,nres
               phi(ii)=phiall(ii,icon)
               theta(ii)=thetall(ii,icon)
@@ -244,6 +262,9 @@ c Produce only a single PDB file for the leading member of the family
               omeg(ii)=omall(ii,icon)
             enddo
             call chainbuild
+
+            endif
+
             cfname=prefixp(:ilen(prefixp))//numm(:ilen(numm))//exten
             OPEN(ipdb,FILE=CFNAME,STATUS='UNKNOWN',FORM='FORMATTED')
 c           print *,'Calling pdbout'