--- /dev/null
+ subroutine read_ref_structure(*)
+C
+C Read the reference structure from the PDB file or from a PDB file or in the form of the dihedral
+C angles.
+C
+ implicit none
+ include 'DIMENSIONS'
+ include 'DIMENSIONS.ZSCOPT'
+ include 'DIMENSIONS.COMPAR'
+ include 'COMMON.IOUNITS'
+ include 'COMMON.GEO'
+ include 'COMMON.VAR'
+ include 'COMMON.INTERACT'
+ include 'COMMON.LOCAL'
+ include 'COMMON.NAMES'
+ include 'COMMON.CHAIN'
+ include 'COMMON.FFIELD'
+ include 'COMMON.SBRIDGE'
+ include 'COMMON.HEADER'
+ include 'COMMON.CONTROL'
+ include 'COMMON.CONTACTS1'
+ include 'COMMON.PEPTCONT'
+ include 'COMMON.TIME1'
+ include 'COMMON.COMPAR'
+ character*4 sequence(maxres)
+ integer rescode
+ double precision x(maxvar)
+ integer itype_pdb(maxres)
+ logical seq_comp
+ integer i,j,k,nres_pdb,iaux
+ double precision ddsc,dist
+ integer ilen
+ external ilen
+C
+ nres0=nres
+ write (iout,*) "pdbref",pdbref
+ if (pdbref) then
+ read(inp,'(a)') pdbfile
+ write (iout,'(2a,1h.)') 'PDB data will be read from file ',
+ & pdbfile(:ilen(pdbfile))
+ open(ipdbin,file=pdbfile,status='old',err=33)
+ goto 34
+ 33 write (iout,'(a)') 'Error opening PDB file.'
+ return1
+ 34 continue
+ do i=1,nres
+ itype_pdb(i)=itype(i)
+ enddo
+ call readpdb(.true.)
+ do i=1,nres
+ iaux=itype_pdb(i)
+ itype_pdb(i)=itype(i)
+ itype(i)=iaux
+ enddo
+ close (ipdbin)
+ nres_pdb=nres
+ nres=nres0
+ nstart_seq=nnt
+ if (nsup.le.(nct-nnt+1)) then
+ do i=0,nct-nnt+1-nsup
+ if (seq_comp(itype(nnt+i),itype_pdb(nstart_sup),
+ & nsup)) then
+ do j=nnt+nsup-1,nnt,-1
+ do k=1,3
+ cref(k,nres+j+i)=cref(k,nres_pdb+j)
+ enddo
+ enddo
+ do j=nnt+nsup-1,nnt,-1
+ do k=1,3
+ cref(k,j+i)=cref(k,j)
+ enddo
+ phi_ref(j+i)=phi_ref(j)
+ theta_ref(j+i)=theta_ref(j)
+ alph_ref(j+i)=alph_ref(j)
+ omeg_ref(j+i)=omeg_ref(j)
+ enddo
+#ifdef DEBUG
+ do j=nnt,nct
+ write (iout,'(i5,3f10.5,5x,3f10.5)')
+ & j,(cref(k,j),k=1,3),(cref(k,j+nres),k=1,3)
+ enddo
+#endif
+ nstart_seq=nnt+i
+ nstart_sup=nnt+i
+ goto 111
+ endif
+ enddo
+ write (iout,'(a)')
+ & 'Error - sequences to be superposed do not match.'
+ return1
+ else
+ do i=0,nsup-(nct-nnt+1)
+ if (seq_comp(itype(nnt),itype_pdb(nstart_sup+i),
+ & nct-nnt+1))
+ & then
+ nstart_sup=nstart_sup+i
+ nsup=nct-nnt+1
+ goto 111
+ endif
+ enddo
+ write (iout,'(a)')
+ & 'Error - sequences to be superposed do not match.'
+ endif
+ 111 continue
+ write (iout,'(a,i5)')
+ & 'Experimental structure begins at residue',nstart_seq
+ else
+ call read_angles(inp,*38)
+ goto 39
+ 38 write (iout,'(a)') 'Error reading reference structure.'
+ return1
+ 39 call chainbuild
+ nstart_sup=nnt
+ nstart_seq=nnt
+ nsup=nct-nnt+1
+ do i=1,2*nres
+ do j=1,3
+ cref(j,i)=c(j,i)
+ enddo
+ enddo
+ endif
+ nend_sup=nstart_sup+nsup-1
+ do i=1,2*nres
+ do j=1,3
+ c(j,i)=cref(j,i)
+ enddo
+ enddo
+ do i=1,nres
+ do j=1,3
+ dc(j,nres+i)=cref(j,nres+i)-cref(j,i)
+ enddo
+ if (itype(i).ne.10) then
+ ddsc = dist(i,nres+i)
+ do j=1,3
+ dc_norm(j,nres+i)=dc(j,nres+i)/ddsc
+ enddo
+ else
+ do j=1,3
+ dc_norm(j,nres+i)=0.0d0
+ enddo
+ endif
+c write (iout,*) "i",i," dc_norm",(dc_norm(k,nres+i),k=1,3),
+c " norm",dc_norm(1,nres+i)**2+dc_norm(2,nres+i)**2+
+c dc_norm(3,nres+i)**2
+ do j=1,3
+ dc(j,i)=c(j,i+1)-c(j,i)
+ enddo
+ ddsc = dist(i,i+1)
+ do j=1,3
+ dc_norm(j,i)=dc(j,i)/ddsc
+ enddo
+ enddo
+c print *,"Calling contact"
+ call contact(.true.,ncont_ref,icont_ref(1,1),
+ & nstart_sup,nend_sup)
+c print *,"Calling elecont"
+ call elecont(.true.,ncont_pept_ref,
+ & icont_pept_ref(1,1),
+ & nstart_sup,nend_sup)
+ write (iout,'(a,i3,a,i3,a,i3,a)')
+ & 'Number of residues to be superposed:',nsup,
+ & ' (from residue',nstart_sup,' to residue',
+ & nend_sup,').'
+ return
+ end