c23456789012345678901234567890123456789012345678901234567890123456789012 program convpdbfi character*3 resnam character*3 restyp(800) character*8 fil(160) character*4 prna character*4 atnam character*6 dummy integer atnum,resnum,resold,maxres,icount(22) integer nres(800),natoms(800),nala(22,800),iala(22,800),iss2(200) integer dist(22,25),resno(160),iss1(200) integer icounf(22,22),totfi(22,22),ibond(22) integer icounk(22),iir(150,800),doublet(22,13,13) integer dt(22,22,18),isumxr(22),isumrx(22),doublet2(22,13,13) integer doublet3(22,13,13),dall(13,13),resno1(800) real x(800,35),y(800,35),z(800,35),totan(22),avgan(22) real mass(800,35),xx(800),yy(800),zz(800) real avgbo(22),ax(6,80,2),ay(6,80,2),az(6,80,2) real elx(800),ely(800),elz(800),elm(800),elmsq(800) real crox(800),croy(800),croz(800),crom(800),cromsq(800) real croxp(800),croyp(800),crozp(800),cromp(800),hyd(21) real theta(150,800),phi(150,800) real bond(22,4700),thet1(22),thet2(22) real sidex(800),sidey(800),sidez(800) real cyss(150,800),blength(150,800) real thetpr(150,800),phip(150,800) c This program calculates virtual bond model coordinates. First c virtual bond lengths are evaluated. Then bond angles will be c determined for each residue. Finally torsional angles between c virtual bonds will be evaluated. We need a single point per c residue, which will be found from the mass center of the points c already considered in evaluating pair radial distribution fns. c These points will be saved as xx, yy and zz, for each residue in c a given protein. Therefore the corresponding array sizes will c be 800. Likewise the alpha carbons will be used for specifying c backbone locations c nala(kk,iiii) designates the number of residues of type kk in c the protein iiii. natoms(kk) designates the total number of c atoms of residue kk, in a given protein, whose coordinates are c listed in PDB. iala(kk,i) denotes the serial index of the ith c residue of type kk along a given protein. fil(1)='1STN.pdb' do 59331 i=1,22 icount(i)=0 ibond(i)=0 do 59332 j=1,22 59332 icounf(i,j)=0 59331 continue c nfile2>1 if more than 1 file is to be read nfile1=1 nfile2=1 c_____________________________________________________________________ do 4293 iiii=nfile1,nfile2 do 40893 j=1,800 xx(j)=0. yy(j)=0. 40893 zz(j)=0. nsb=0 conv=180./3.141592654 write(6,*)'********',iiii open(10,file=fil(iiii)) do 3 i=1,600 do 3 j=1,35 x(i,j)=0. y(i,j)=0. 3 z(i,j)=0. i=0 resold=1 5 read(10,50)dummy 50 Format(A6) 40 format (A4) 55 FORMAT (A6,1X,I4,1X,A4,1X,A3,1X,A1,1X,I3,5X,3F8.3,5X,F6.2) 54 format (A6,1X,I3,1X,A1,2X,I3,2X,13(A3,1X),1X,A4) 56 format (A6,1X,I3,1X,A3,1X,A1,2X,I3,4X,A3,1X,A1,2X,I3) if (dummy.eq.'SEQRES') then backspace (10) read(10,54)dummy,nrow,s,maxres,(restyp(13*(nrow-1)+j),j=1,13),prna write(70,54)dummy,nrow,s,maxres,(restyp(13*(nrow-1)+j),j=1,13), + prna endif maxres=maxresold if(dummy.ne.'ATOM ') goto 5 backspace (10) 12 i=i+1 1231 read(10,50)dummy if(dummy.eq.'HETATM') goto 33 if(dummy.eq.'SIGATM') goto 1231 backspace (10) read(10,55)dummy,atnum,atnam,resnam,s,resnum,x(resnum,i), + y(resnum,i),z(resnum,i),bf write(70,55)dummy,atnum,atnam,resnam,s,resnum, + x(resnum,i),y(resnum,i),z(resnum,i),bf if(atnam.eq.' CA ') then write(92,*)resnum,bf endif if(atnum.eq.1) resno1(iiii)=resnum if(resnum.ne.resold)then natoms(resold)=i-1 x(resnum,1)=x(resnum,i) y(resnum,1)=y(resnum,i) z(resnum,1)=z(resnum,i) restyp(resnum)=resnam i=1 resold=resnum endif if(dummy.eq.'TER ') goto 33 c if(atnum.eq.757) goto 33 go to 12 c************************************************************* 33 maxat=atnum c write(91,*)'resnum',resnum natoms(resnum)=i-1 resno(iiii)=resnum do 29811 i=1,resnum blength(iiii,i)=0. theta(iiii,i)=0. thetpr(iiii,i)=0. phip(iiii,i)=0. 29811 phi(iiii,i)=0. kala=1 kgly=1 kval=1 kleu=1 kile=1 kphe=1 karg=1 klys=1 kglu=1 kasp=1 kasn=1 kgln=1 kpro=1 khis=1 ktyr=1 ktrp=1 kthr=1 kser=1 kcys=1 kmet=1 kcyss=1 do 19001 i=resno1(iiii),resnum write(91,*)i,restyp(i) c atomic mass of different atoms N, O, C.. mass(i,1)=15. mass(i,2)=13. mass(i,3)=12. mass(i,4)=16. do 12022 j=5,14 12022 mass(i,j)=12. do 12029 j=15,25 12029 mass(i,j)=1. if (restyp(i).eq.'ALA') THEN iir(iiii,i)=2 if (natoms(i).ge.5) then write(91,*)'kkkkkkk',i iala(2,kala)=i kala=kala+1 xx(i)=x(i,5) yy(i)=y(i,5) zz(i)=z(i,5) endif endif if (restyp(i).eq.'VAL') THEN iir(iiii,i)=3 if (natoms(i).ge.7) then write(91,*)'kkkkkkk',i iala(3,kval)=i kval=kval+1 xx(i)=(x(i,6)+x(i,7))/2.0 yy(i)=(y(i,6)+y(i,7))/2.0 zz(i)=(z(i,6)+z(i,7))/2.0 endif endif if (restyp(i).eq.'LEU') THEN iir(iiii,i)=5 if (natoms(i).ge.8) then write(91,*)'kkkkkkk',i iala(5,kleu)=i kleu=kleu+1 xx(i)=(x(i,8)+x(i,7))/2.0 yy(i)=(y(i,8)+y(i,7))/2.0 zz(i)=(z(i,8)+z(i,7))/2.0 endif endif if (restyp(i).eq.'ILE') THEN iir(iiii,i)=4 if (natoms(i).ge.8) then write(91,*)'kkkkkkk',i iala(4,kile)=i kile=kile+1 xx(i)=x(i,8) yy(i)=y(i,8) zz(i)=z(i,8) endif endif if (restyp(i).eq.'PHE') THEN iir(iiii,i)=16 mass(i,2)=17. if (natoms(i).ge.11) then write(91,*)'kkkkkkk',i iala(16,kphe)=i kphe=kphe+1 do 48331 iji=6,11 xx(i)=xx(i)+x(i,iji) yy(i)=yy(i)+y(i,iji) 48331 zz(i)=zz(i)+z(i,iji) xx(i)=xx(i)/6. yy(i)=yy(i)/6. zz(i)=zz(i)/6. endif endif if (restyp(i).eq.'CYS') THEN do 38999 ij=1,nsb if(i.eq.iss1(ij).or.i.eq.iss2(ij))then iala(21,kcyss)=i cyss(iiii,i)=1. kcyss=kcyss+1 endif 38999 continue iala(14,kcys)=i kcys=kcys+1 mass(i,6)=33. iir(iiii,i)=14 if (natoms(i).ge.6) then write(91,*)'kkkkkkk',i yy(i)=y(i,6) zz(i)=z(i,6) xx(i)=x(i,6) endif endif if (restyp(i).eq.'MET') THEN iir(iiii,i)=15 mass(i,7)=32. if (natoms(i).ge.7) then write(91,*)'kkkkkkk',i iala(15,kmet)=i kmet=kmet+1 xx(i)=x(i,7) yy(i)=y(i,7) zz(i)=z(i,7) endif endif if (restyp(i).eq.'HIS') then iir(iiii,i)=19 mass(i,7)=14. mass(i,10)=14. if (natoms(i).ge.10) then write(91,*)'kkkkkkk',i iala(19,khis)=i khis=khis+1 do 48332 iji=6,10 xx(i)=xx(i)+x(i,iji) yy(i)=yy(i)+y(i,iji) 48332 zz(i)=zz(i)+z(i,iji) xx(i)=xx(i)/5. yy(i)=yy(i)/5. zz(i)=zz(i)/5. endif endif if (restyp(i).eq.'TRP') then iir(iiii,i)=18 mass(i,9)=14. if (natoms(i).ge.14) then write(91,*)'kkkkkkk',i iala(18,ktrp)=i ktrp=ktrp+1 do 48333 iji=7,14 xx(i)=xx(i)+x(i,iji) yy(i)=yy(i)+y(i,iji) 48333 zz(i)=zz(i)+z(i,iji) xx(i)=xx(i)/8. yy(i)=yy(i)/8. zz(i)=zz(i)/8. endif endif if (restyp(i).eq.'ASP') then iir(iiii,i)=8 mass(i,7)=16. mass(i,8)=16. if (natoms(i).ge.8) then write(91,*)'kkkkkkk',i iala(8,kasp)=i kasp=kasp+1 xx(i)=(x(i,8)+x(i,7))/2.0 yy(i)=(y(i,8)+y(i,7))/2.0 zz(i)=(z(i,8)+z(i,7))/2.0 endif endif if (restyp(i).eq.'GLU') then iir(iiii,i)=10 mass(i,9)=16. mass(i,8)=16. if (natoms(i).ge.9) then write(91,*)'kkkkkkk',i iala(10,kglu)=i kglu=kglu+1 xx(i)=(x(i,8)+x(i,9))/2.0 yy(i)=(y(i,8)+y(i,9))/2.0 zz(i)=(z(i,8)+z(i,9))/2.0 endif endif if (restyp(i).eq.'ASN') then iir(iiii,i)=9 mass(i,7)=16. mass(i,8)=14. if (natoms(i).ge.8) then write(91,*)'kkkkkkk',i iala(9,kasn)=i kasn=kasn+1 xx(i)=(x(i,8)+x(i,7))/2.0 yy(i)=(y(i,8)+y(i,7))/2.0 zz(i)=(z(i,8)+z(i,7))/2.0 endif endif if (restyp(i).eq.'GLN') then iir(iiii,i)=11 mass(i,8)=16. mass(i,9)=14. if (natoms(i).ge.9) then write(91,*)'kkkkkkk',i iala(11,kgln)=i kgln=kgln+1 xx(i)=(x(i,8)+x(i,9))/2.0 yy(i)=(y(i,8)+y(i,9))/2.0 zz(i)=(z(i,8)+z(i,9))/2.0 endif endif if (restyp(i).eq.'SER') then iir(iiii,i)=6 mass(i,6)=16. if (natoms(i).ge.6) then write(91,*)'kkkkkkk',i iala(6,kser)=i kser=kser+1 xx(i)=x(i,6) yy(i)=y(i,6) zz(i)=z(i,6) endif endif if (restyp(i).eq.'THR') then iir(iiii,i)=7 mass(i,6)=16. if (natoms(i).ge.6) then write(91,*)'kkkkkkk',i iala(7,kthr)=i kthr=kthr+1 xx(i)=x(i,6) yy(i)=y(i,6) zz(i)=z(i,6) endif endif if (restyp(i).eq.'LYS') then iir(iiii,i)=12 mass(i,9)=15. if (natoms(i).ge.9) then write(91,*)'kkkkkkk',i iala(12,klys)=i klys=klys+1 xx(i)=x(i,9) yy(i)=y(i,9) zz(i)=z(i,9) endif endif if (restyp(i).eq.'ARG') then iir(iiii,i)=13 mass(i,8)=15. mass(i,10)=16. mass(i,11)=16. if (natoms(i).ge.11) then write(91,*)'kkkkkkk',i iala(13,karg)=i karg=karg+1 xx(i)=(x(i,8)+x(i,10)+x(i,11))/3.0 yy(i)=(y(i,8)+y(i,10)+y(i,11))/3.0 zz(i)=(z(i,8)+z(i,10)+z(i,11))/3.0 endif endif if (restyp(i).eq.'TYR') then iir(iiii,i)=17 mass(i,12)=17. if (natoms(i).ge.12) then write(91,*)'kkkkkkk',i iala(17,ktyr)=i ktyr=ktyr+1 do 48334 iji=6,12 xx(i)=xx(i)+x(i,iji) yy(i)=yy(i)+y(i,iji) 48334 zz(i)=zz(i)+z(i,iji) xx(i)=xx(i)/7. yy(i)=yy(i)/7. zz(i)=zz(i)/7. endif endif if (restyp(i).eq.'GLY') then iir(iiii,i)=1 if (natoms(i).ge.1) then write(91,*)'kkkkkkk',i iala(1,kgly)=i kgly=kgly+1 xx(i)=x(i,2) yy(i)=y(i,2) zz(i)=z(i,2) endif endif if (restyp(i).eq.'PRO') then iir(iiii,i)=20 if (natoms(i).ge.7) then write(91,*)'kkkkkkk',i iala(20,kpro)=i kpro=kpro+1 xx(i)=(x(i,5)+x(i,6)+x(i,7))/3.0 yy(i)=(y(i,5)+y(i,6)+y(i,7))/3.0 zz(i)=(z(i,5)+z(i,6)+z(i,7))/3.0 endif endif 19001 continue sum=0. sum=kala+kgly+kval+kleu+kile+kphe+karg+klys+kglu+kasp+ + kasn+kgln+kpro+khis+ktyr+ktrp+kthr+kser+kcys+kmet c write(91,*)'sum',sum kala=kala-1 kgly=kgly-1 kval=kval-1 kleu=kleu-1 kile=kile-1 kphe=kphe-1 karg=karg-1 klys=klys-1 kglu=kglu-1 kasp=kasp-1 kasn=kasn-1 kgln=kgln-1 kpro=kpro-1 khis=khis-1 ktyr=ktyr-1 ktrp=ktrp-1 kthr=kthr-1 kser=kser-1 kcys=kcys-1 kmet=kmet-1 kcyss=kcyss-1 nala(2,iiii)=kala nala(1,iiii)=kgly nala(3,iiii)=kval nala(4,iiii)=kile nala(5,iiii)=kleu nala(6,iiii)=kser nala(7,iiii)=kthr nala(8,iiii)=kasp nala(9,iiii)=kasn nala(10,iiii)=kglu nala(11,iiii)=kgln nala(12,iiii)=klys nala(13,iiii)=karg nala(14,iiii)=kcys nala(15,iiii)=kmet nala(16,iiii)=kphe nala(17,iiii)=ktyr nala(18,iiii)=ktrp nala(19,iiii)=khis nala(20,iiii)=kpro nala(21,iiii)=kcyss sum=0 do 89 k=1,20 sum=sum+nala(k,1) 89 continue write(*,*)'number of sum',sum,'unrecorded',maxres-sum c__________________________________________________________________ c Residues with unknown coords are excluded from calculations c__________________________________________________________________ c CALCULATION OF BACKBONE BOND LENGTHS AND ANGLES DO 50572 K=resno1(iiii)+1,resnum CROX(K) = 0. CROY(K) = 0. CROZ(K) = 0. ELX(K) = 0. ELY(K) = 0. 50572 ELZ(K) = 0. DO 50571 K=resno1(iiii)+1,resnum if (x(K,2).eq.0.or.x(K-1,2).eq.0)then write(91,*)'kkkkkkkkk',k,iir(iiii,k) goto 50571 endif ELX(K) = x(K,2) - x(K-1,2) ELY(K) = y(K,2) - y(K-1,2) ELZ(K) = z(K,2) - z(K-1,2) ELMSQ(K) = ELX(K)*ELX(K) + ELY(K)*ELY(K) + ELZ(K)*ELZ(K) write(91,*)sqrt(ELMSQ(K)),k 50571 ELM(K) = SQRT (ELMSQ(K)) c write(6,*)iiii,elm(10),elm(20),elm(30),elm(40) if (x(3,2).eq.0.or.x(2,2).eq.0) goto 50576 ELX(1) = x(3,2) - x(2,2) ELY(1) = y(3,2) - y(2,2) ELZ(1) = z(3,2) - z(2,2) ELMSQ(1) = ELX(1)*ELX(1) + ELY(1)*ELY(1) + ELZ(1)*ELZ(1) ELM(1) = SQRT (ELMSQ(1)) 50576 k=resnum-1 kf=resnum+1 if (x(k,2).eq.0.or.x(k-1,2).eq.0) goto 50577 ELX(kf) = x(k,2) - x(k-1,2) ELY(kf) = y(k,2) - y(k-1,2) ELZ(kf) = z(k,2) - z(k-1,2) ELMSQ(kf) = ELX(kf)*ELX(kf) + ELY(kf)*ELY(kf) + ELZ(kf)*ELZ(kf) ELM(kf) = SQRT (ELMSQ(kf)) 50577 continue c__________________________________________________________________ c write(2,*)'CALCULATION OF VIRTUAL BOND ANGLES' c WRITE(2,*)' ' c here theta(iiii,mi) means the bond angle of the mi th residue c in the iiii th protein.phi(iiii,k) c has similar meaning as theta(iiii,mi). do 10000 jres1=1,20 DO 40571 i=1,nala(jres1,iiii) mi=iala(jres1,i) if (mi.lt.2) goto 40571 if (mi.gt.(resnum-1)) goto 40571 if (elx(mi).eq.0.or.elx(mi+1).eq.0) then write(6,*)'one zero',iiii,mi,jres1 goto 40571 endif icount(jres1)=icount(jres1)+1 kc=icount(jres1) dot=elx(mi+1)*elx(mi)+ely(mi+1)*ely(mi) dot = -(dot + elz(mi+1)*elz(mi))/elm(mi)/elm(mi+1) theta(iiii,mi)=acos(dot)*conv write(91,*)'th',theta(iiii,mi),mi 40571 continue totan(jres1)=icount(jres1) 10000 continue c__________________________________________________________________ C **CROSS PRODUCTS OF BOND VECTORS I AND (I+1) AND MAGNITUDES** c here we include the hypothetical first and last bonds as well. DO 801 K = resno1(iiii),resnum if (elx(k).ne.0.and.elx(k+1).ne.0) then croX(K) = ELY(K)*ELZ(K+1)-ELZ(K)*ELY(K+1) croY(K) = -ELX(K)*ELZ(K+1)+ELZ(K)*ELX(K+1) croZ(K) = ELX(K)*ELY(K+1)-ELY(K)*ELX(K+1) croM(K) = crox(K)*crox(k)+croy(K)*croy(k)+croz(K)*croz(k) crom(k)=crom(k)**0.5 c write(6,*)elx(k),ely(k),elz(k),elx(k+1),ely(k+1),elz(k+1) c write(6,*)'cros pr',iiii,k,crox(k),croy(k),croz(k),crom(k) endif 801 CONTINUE c write(2,*)'CALCULATION OF TORSIONAL ANGLES' DO 804 K = resno1(iiii),resnum if (crox(k).ne.0.and.crox(k-1).ne.0) then ara=CROX(K)*CROX(K-1)+CROY(K)*CRoy(k-1)+CROZ(K)*CROZ(K-1) EKFY = -1./CROM(K)/CROM(K-1) cosfi = ara*EKFY if (cosfi.gt.1.) cosfi=1.0 if (cosfi.lt.-1.) cosfi=-1.0 FI = ACOS(cosfi) gpm = crox(k-1)*elx(k+1)+croy(k-1)*ely(k+1)+croz(k-1)*elz(k+1) if (gpm.gt.0.0) fi=-fi faydeg = fi*conv c write(6,*)'angle',iiii,k,faydeg,cosfi,ara i1=iir(iiii,k-1) i2=iir(iiii,k) icounf(i1,i2)=icounf(i1,i2)+1 kc=icounf(i1,i2) phi(iiii,k)=faydeg c write(6,*)' ' write(91,*)'fay',phi(iiii,k),k endif phi(iiii,2)=0. 804 CONTINUE c__________________________________________________________________ c Now we know each bond bending, torsion, location of side groups c in space. Locations wrt backbone will be described by three c variables: these are one length and two orientations as follows: c (1)length of side group (i.e.distance between alpha-C c and center xx,yy,zz of side group), (2) The angle the side group c i makes with the plane of the alpha carbons i-1, i and i+1, c and (3) The angle the side group makes with the normal plane to c bisecting the bond angle i. do 10032 k=resno1(iiii),resnum if (restyp(k).eq.'GLY') then sidex(k)=0. sidey(k)=0. sidez(k)=0. goto 10032 endif if (xx(k).ne.0.0) then x1=xx(k)-x(k,2) y1=yy(k)-y(k,2) z1=zz(k)-z(k,2) ara = x1*x1 + y1*y1 + z1*z1 ara=SQRT(ara) blength(iiii,k)=ara jres=iir(iiii,k) ibond(jres)=ibond(jres)+1 mco=ibond(jres) bond(jres,mco)=ara sidex(k)=x1/bond(jres,mco) sidey(k)=y1/bond(jres,mco) sidez(k)=z1/bond(jres,mco) c write(6,*)iiii,k,jres,ibond(jres),bond(jres,mco) endif 10032 continue do 10001 jres1=1,20 DO 42571 i=1,nala(jres1,iiii) mi=iala(jres1,i) if (restyp(mi).eq.'GLY') THEN thetpr(iiii,mi)=0. phip(iiii,mi)=0. goto 42571 endif if (elx(mi).eq.0.or.sidex(mi).eq.0) then write(6,*)'one side group zero',iiii,mi,jres1 goto 42571 endif dot=sidex(mi)*elx(mi)+sidey(mi)*ely(mi) dot = -(dot + sidez(mi)*elz(mi))/elm(mi) thetpr(iiii,mi)=acos(dot)*conv 42571 continue 10001 continue if (elx(1).eq.0.or.sidex(1).eq.0) then write(6,*)'one terminal side group zero' goto 42574 endif dot=-sidex(1)*elx(2)-sidey(1)*ely(2) dot = -(dot - sidez(1)*elz(2))/elm(2) thetpr(iiii,1)=acos(dot)*conv 42574 continue c cross product of bond i and side bond appended to residue i DO 80117 K = resno1(iiii)+1,resnum if (elx(k).ne.0.and.sidex(k).ne.0) then croxp(K) = ELY(K)*sideZ(K)-ELZ(K)*sideY(K) c write(6,*)'now',croxp(k),ely(k),elz(k),sidez(k),sidey(k) croyp(K) = -ELX(K)*sideZ(K)+ELZ(K)*sideX(K) croZp(K) = ELX(K)*sideY(K)-ELY(K)*sideX(K) croMp(K) = croxp(K)*croxp(k)+croyp(K)*croyp(k)+crozp(K)*crozp(k) cromp(k)=cromp(k)**0.5 c write(6,*)'cros pr',iiii,k,croxp(k),croyp(k),crozp(k),cromp(k) endif 80117 CONTINUE if (elx(2).ne.0.and.sidex(1).ne.0) then croxp(1) = -ELY(2)*sideZ(1)+ELZ(2)*sideY(1) croyp(1) = +ELX(2)*sideZ(1)-ELZ(2)*sideX(1) croZp(1) = -ELX(2)*sideY(1)+ELY(2)*sideX(1) croMp(1) = croxp(1)*croxp(1)+croyp(1)*croyp(1)+crozp(1)*crozp(1) cromp(1)=cromp(1)**0.5 endif c write(2,*)'CALCULATION OF SIDEGROUP TORSIONAL ANGLES' do 8104 k=resno1(iiii)+1,resnum if (restyp(i).eq.'GLY') goto 8104 if (croxp(k).ne.0.and.crox(k-1).ne.0) then ara=CROXp(K)*CROX(K-1)+CROYp(K)*CRoy(k-1)+CROZp(K)*CROZ(K-1) EKFY = -1./CROMp(K)/CROM(K-1) cosfi = ara*EKFY if (cosfi.gt.1.) cosfi=1.0 if (cosfi.lt.-1.) cosfi=-1.0 FI = ACOS(cosfi) gpm = crox(k-1)*sidex(k)+croy(k-1)*sidey(k)+croz(k-1)*sidez(k) if (gpm.gt.0.0) fi=-fi faydeg = fi*conv phip(iiii,k)=faydeg if (k.eq.10) then write(6,*)'bond k-1=',elx(k-1),ely(k-1),elz(k-1) write(6,*)'bond k=',elx(k),ely(k),elz(k) write(6,*)'bond k+1=',elx(k+1),ely(k+1),elz(k+1) write(6,*)'side k=',sidex(k),sidey(k),sidez(k) write(6,*)'crox(k-1)',crox(k-1),croy(k-1),croz(k-1) write(6,*)'crox(k)',crox(k),croy(k),croz(k) write(6,*)'croxp(k)',croxp(k),croyp(k),crozp(k) write(6,*)'bond k+1=',elx(k+1),ely(k+1),elz(k+1) write(6,*)'torsions:',phi(iiii,k),phip(iiii,k) write(6,*)'bond angles:',theta(iiii,k),thetpr(iiii,k) endif endif 8104 CONTINUE if (restyp(1).eq.'GLY') goto 8105 if (croxp(1).ne.0.and.crox(2).ne.0) then ara=CROXp(1)*CROX(2)+CROYp(1)*CRoy(2)+CROZp(1)*CROZ(2) EKFY = 1./CROMp(1)/CROM(2) cosfi = ara*EKFY if (cosfi.gt.1.) cosfi=1.0 if (cosfi.lt.-1.) cosfi=-1.0 FI = ACOS(cosfi) gpm = crox(2)*sidex(1)+croy(2)*sidey(1)+croz(2)*sidez(1) if (gpm.gt.0.0) fi=-fi faydeg = fi*conv phip(iiii,1)=faydeg endif 8105 continue c We also considered terminal side groups. For that purpose we c assumed that there are two hypothetical backbone bonds, one in c the same direction as the second virtual bond (elx(3),ely(3),etc) c such that the first virtual bond has a virtual torsion of 0 deg, c and the other same direction as bond n-1, such that the terminal c bond is trans. These are called elx-y-z(1) and elx-y-z(kf), c kf=resnum+1, and defined in the section following loop 50571. c*********************************************************************** if(iiii.eq.2) then do 48933 i=1,resno(iiii) ax(iiii,i,1)=x(i,2) ay(iiii,i,1)=y(i,2) az(iiii,i,1)=z(i,2) ax(iiii,i,2)=xx(i) ay(iiii,i,2)=yy(i) 48933 az(iiii,i,2)=zz(i) endif open(7,file='cor1cho1.out1') C COLUMNS: 1 = residue type c 2-4: x-y-z coordinates of alpha carbons c 5-7: x-y-z coordinates of side group centroids write(7,*)resno(iiii)-5 c write(7,*)'backbone and side group coordinates of ',fil(m) do 20000 i=6,resno(iiii) 20000 write(7,78256)iir(iiii,i),x(i,2),y(i,2),z(i,2), :xx(i),yy(i),zz(i) 4293 continue close(7) c______________________________________________________________________ open(7,file='fi1cho1.out1') C COLUMNS: 1=residue type c 2=torsional angle of preceding bond, 3=bond angle, c 4=side group torsion expressed by preceding bond, c 5=side group bond length, 6=side group bond angle do 27990 m=nfile1,nfile2 write(7,*)resno(m) do 27000 i=7,resno(m) if(cyss(m,i).eq.1.) iir(m,i)=21 27000 write(7,78255)iir(m,i),phi(m,i),elm(i),theta(m,i) :,phip(m,i), :blength(m,i),thetpr(m,i) 27990 continue close(7) 78255 format(i5,6f9.2) 78256 format(i5,6f9.2) 7824 format(3i5,4f13.4) c______________________________________________________________________ stop end