Merge branch 'lipid' of mmka.chem.univ.gda.pl:unres into lipid
[unres.git] / source / unres / src_MD-M-newcorr / rmsd.F
diff --git a/source/unres/src_MD-M-newcorr/rmsd.F b/source/unres/src_MD-M-newcorr/rmsd.F
new file mode 100644 (file)
index 0000000..82749b4
--- /dev/null
@@ -0,0 +1,172 @@
+      subroutine rms_nac_nnc(rms,frac,frac_nn,co,lprn)
+        implicit real*8 (a-h,o-z)
+        include 'DIMENSIONS'
+        include 'COMMON.CHAIN'
+        include 'COMMON.CONTACTS'
+        include 'COMMON.IOUNITS'
+        double precision przes(3),obr(3,3)
+        logical non_conv,lprn
+c        call fitsq(rms,c(1,nstart_seq),cref(1,nstart_sup),nsup,przes,
+c     &             obr,non_conv)
+c        rms=dsqrt(rms)
+        call rmsd(rms)
+        call contact(.false.,ncont,icont,co)
+        frac=contact_fract(ncont,ncont_ref,icont,icont_ref)
+        frac_nn=contact_fract_nn(ncont,ncont_ref,icont,icont_ref)
+        if (lprn) write (iout,'(a,f8.3/a,f8.3/a,f8.3/a,f8.3)')
+     &    'RMS deviation from the reference structure:',rms,
+     &    ' % of native contacts:',frac*100,
+     &    ' % of nonnative contacts:',frac_nn*100,
+     &    ' contact order:',co
+
+      return
+      end      
+c---------------------------------------------------------------------------
+      subroutine rmsd(drms)
+      implicit real*8 (a-h,o-z)
+      include 'DIMENSIONS'
+#ifdef MPI
+      include 'mpif.h'
+#endif
+      include 'COMMON.CHAIN'
+      include 'COMMON.IOUNITS'  
+      include 'COMMON.INTERACT'
+      include 'COMMON.CONTROL'
+      logical non_conv
+      double precision przes(3),obrot(3,3)
+      double precision ccopy(3,maxres2+2),crefcopy(3,maxres2+2)
+      nperm=1
+      do i=1,symetr
+      nperm=nperm*i
+      enddo
+      iatom=0
+      rminroz=100d2
+c      print *,"nz_start",nz_start," nz_end",nz_end
+c      if (symetr.le.1) then
+      do kkk=1,nperm
+c      do i=nz_start,nz_end
+c        iatom=iatom+1
+c        iti=itype(i)
+c        do k=1,3
+c         ccopy(k,iatom)=c(k,i+nstart_seq-nstart_sup)
+c         crefcopy(k,iatom,kkk)=cref(k,i,kkk)
+c        enddo
+c        if (iz_sc.eq.1.and.iti.ne.10) then
+c          iatom=iatom+1
+c          do k=1,3
+c           ccopy(k,iatom)=c(k,nres+i+nstart_seq-nstart_sup)
+c           crefcopy(k,iatom,kkk)=cref(k,nres+i,kkk)
+c          enddo
+c        endif
+c      enddo
+c      else
+c      do kkk=1,nperm
+      iatom=0
+      do i=nz_start,nz_end
+        iatom=iatom+1
+        iti=itype(i)
+        do k=1,3
+         ccopy(k,iatom)=c(k,i+nstart_seq-nstart_sup)
+         crefcopy(k,iatom)=cref(k,i,kkk)
+        enddo
+        if (iz_sc.eq.1.and.iti.ne.10) then
+          iatom=iatom+1
+          do k=1,3
+           ccopy(k,iatom)=c(k,nres+i+nstart_seq-nstart_sup)
+           crefcopy(k,iatom)=cref(k,nres+i,kkk)
+          enddo
+        endif
+      enddo
+c      enddo
+c      endif
+      
+c ----- diagnostics
+c         do kkk=1,nperm
+c          write (iout,*) 'Ccopy and CREFcopy'
+c          print '(i5,3f10.5,5x,3f10.5)',(k,(ccopy(j,k),j=1,3),
+c     &           (crefcopy(j,k),j=1,3),k=1,iatom)
+c          write (iout,'(i5,3f10.5,5x,3f10.5)') (k,(ccopy(j,k),j=1,3),
+c     &           (crefcopy(j,k),j=1,3),k=1,iatom)
+c         enddo
+c ----- end diagnostics
+c      do kkk=1,nperm
+      call fitsq(roznica,ccopy(1,1),crefcopy(1,1),iatom,
+     &                                      przes,obrot,non_conv) 
+      if (non_conv) then
+          print *,'Problems in FITSQ!!! rmsd'
+          write (iout,*) 'Problems in FITSQ!!! rmsd'
+          print *,'Ccopy and CREFcopy'
+          write (iout,*) 'Ccopy and CREFcopy'
+          print '(i5,3f10.5,5x,3f10.5)',(k,(ccopy(j,k),j=1,3),
+     &           (crefcopy(j,k),j=1,3),k=1,iatom)
+          write (iout,'(i5,3f10.5,5x,3f10.5)') (k,(ccopy(j,k),j=1,3),
+     &           (crefcopy(j,k),j=1,3),k=1,iatom)
+#ifdef MPI
+c          call mpi_abort(mpi_comm_world,ierror,ierrcode)
+           roznica=100.0d10
+#else          
+          stop
+#endif
+       endif
+c       write (iout,*) "roznica", roznica,kkk
+       if (roznica.le.rminroz) rminroz=roznica
+       enddo
+       drms=dsqrt(dabs(rminroz))
+c ---- diagnostics
+c        write (iout,*) "nperm,symetr", nperm,symetr
+c ---- end diagnostics
+       return
+       end
+
+c--------------------------------------------
+      subroutine rmsd_csa(drms)
+      implicit real*8 (a-h,o-z)
+      include 'DIMENSIONS'
+#ifdef MPI
+      include 'mpif.h'
+#endif
+      include 'COMMON.CHAIN'
+      include 'COMMON.IOUNITS'  
+      include 'COMMON.INTERACT'
+      logical non_conv
+      double precision przes(3),obrot(3,3)
+      double precision ccopy(3,maxres2+2),crefcopy(3,maxres2+2)
+      kkk=1
+      iatom=0
+      do i=nz_start,nz_end
+        iatom=iatom+1
+        iti=itype(i)
+        do k=1,3
+         ccopy(k,iatom)=c(k,i)
+         crefcopy(k,iatom)=crefjlee(k,i)
+        enddo
+        if (iz_sc.eq.1.and.iti.ne.10) then
+          iatom=iatom+1
+          do k=1,3
+           ccopy(k,iatom)=c(k,nres+i)
+           crefcopy(k,iatom)=crefjlee(k,nres+i)
+          enddo
+        endif
+      enddo
+
+      call fitsq(roznica,ccopy(1,1),crefcopy(1,1),iatom,
+     &                                      przes,obrot,non_conv) 
+      if (non_conv) then
+          print *,'Problems in FITSQ!!! rmsd_csa'
+          write (iout,*) 'Problems in FITSQ!!! rmsd_csa'
+          print *,'Ccopy and CREFcopy'
+          write (iout,*) 'Ccopy and CREFcopy'
+          print '(i5,3f10.5,5x,3f10.5)',(k,(ccopy(j,k),j=1,3),
+     &           (crefcopy(j,k),j=1,3),k=1,iatom)
+          write (iout,'(i5,3f10.5,5x,3f10.5)') (k,(ccopy(j,k),j=1,3),
+     &           (crefcopy(j,k),j=1,3),k=1,iatom)
+#ifdef MPI
+          call mpi_abort(mpi_comm_world,ierror,ierrcode)
+#else          
+          stop
+#endif
+       endif
+       drms=dsqrt(dabs(roznica))
+       return
+       end
+