C FREEHEL98.FOR 1 June 1998 C This is the Fortran program code for helix analysis program FREEHELIX. C For operating instructions, see the separate file FREEHEL98.TEX. This C program is radically changed from the earlier NEWHELIX, deleting the C old RADJ and TADJ, and adding the helix-independent parameters VALL C (total bending), VROL (roll), VTIL (tilt), VTWI (twist), VSLI (slide), C VRIS (rise) and VSHF (shift). This program can be used with DNA helices C of any degree of bending. C Richard E. Dickerson, UCLA C ********************************************************************* C THE FREEHELIX PROGRAM BEGINS HERE: C ********************************************************************* PROGRAM FREEHELIX CHARACTER * 4 NAME,MH CHARACTER * 1 IR,TITL CHARACTER * 60 FMT COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),ICORL X ,L15 COMMON/INPUT/CELL(6),FMT,NPMIN,NPMAX,IFPUN,IFLAGP X ,IHELIX,IBROLL,ICYLIN,ITORNG COMMON/ATOM/NATM,ATM(3,2000),NAME(2,2000) COMMON/BROOK/IBARF,IFRES,NBP CALL INIT CALL XDRED IF (IHELIX.NE.0) CALL HELIX IF (IBROLL.NE.0) CALL BROLL IF (ICYLIN.NE.0) CALL CYLIN IF (ITORNG.NE.0) CALL TORANG C RETURN END C------------------------- C SUBROUTINE ANGLE C------------------------- C SUGAR BOND ANGLES CALCULATING SUBROUTINE SUBROUTINE ANGLE (NA,TE,ATM) DIMENSION ATM(3,2000),NA(3) DIMENSION AV(3),BV(3),CV(3),AA(3) DO 10 J=1,3 J1=J+1 IF (J1.GT.3) J1=1 NS=NA(J) NE=NA(J1) AV(J)=ATM(1,NE)-ATM(1,NS) BV(J)=ATM(2,NE)-ATM(2,NS) CV(J)=ATM(3,NE)-ATM(3,NS) 10 CONTINUE DO 20 I=1,3 AA(I)=AV(I)**2+BV(I)**2+CV(I)**2 20 CONTINUE TE= ACOS((AA(1)+AA(2)-AA(3))/(2.*SQRT(AA(1)*AA(2)))) RETURN END C------------------------- C SUBROUTINE BROLL C------------------------- SUBROUTINE BROLL CHARACTER * 4 MH CHARACTER * 1 IR,TITL,NTP,NAT CHARACTER * 4 ILAB,XPS,XLS,XPL INTEGER * 2 IRES DOUBLE PRECISION CO1,SI1,CO2,SI2,A1R,B1R,A2R,B2R COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),ICORL COMMON/ATOM/NATM,ATM(3,2000),NTP(2000),IRES(2000),NAT(5,2000) COMMON/BROLLP/INPBRL,NC(2,50),NB1(2,50),NB2(2,50),COSANG(50,50) COMMON/BROOK/IBARF,IFRES,NPLA DIMENSION A(3,3),EIG(20,3),NA1(15),PHI(50,3,2),THE(50,3,4),TIL(50) * ,NA2(15),XMIN(30,3),XMAJ(30,3),PLA(50,3,3),ANG(3),SV(3),RV(2,3) * ,PV(2,3),PM(2),PTW(50),ALB(2),BKL(50),SLI(50),DIS(50),HYP(50) * ,ROT12(50),AM(3),BM(3),SLJ(50),SIG(12),SSQ(12),SD(12),AV(12) * ,TV(2,3),TV2(2),ZALL(50),ZTIL(50),ZROL(50),ZSLI(50),ZTWI(50) * ,CUP(50),ZRIS(50),ZSLP(50),IANGALL(50,50),NJAY(50) 500 FORMAT(1X,' ROLL+TILT OUTPUT, FREEHELIX98') 501 FORMAT(/,I5,' ATOMS AND ',I5,' BASE PAIRS',I5) 502 FORMAT (/,' ',120A1) 503 FORMAT(//,' BASE PAIR',I5,/) 504 FORMAT(2I5) 505 FORMAT(15I5) 506 FORMAT(1H1) 507 FORMAT(//,' STRAND',I4,' BASE NORMAL COSINES AND ANGLES',/) 508 FORMAT(3F10.5,5X,3F10.2) 509 FORMAT(//,' STRAND 1 ROLL AND TILT ANGLES STRAND 2 ROLL A *ND TILT ANGLES') C........1.........2.........3.........4.........5.........6.........712 510 FORMAT(' TIP INCL ROLL TILT RADJ TADJ TIP INCL ROLL * TILT RADJ TADJ') 511 FORMAT(//,' BEST PLANE THROUGH BOTH BASES') C........1.........2.........3.........4.........5.........6.........712 512 FORMAT(' VALL VTIL VROL VSLI VTWI VRIS VSHF INC *L TIP TILT ROLL SLIDE CUP PROP BUCK X DSP Y DSP') 4513 FORMAT(1X,6F6.2,2X,6F6.2) 513 FORMAT(1X,7F7.2,2X,10F7.2) 514 FORMAT(11F7.2) 516 FORMAT(1H ) 517 FORMAT(' NOTE: Angles are calculated from 5" end to 3" end of str *and 1, and signs') 518 FORMAT(' of angles also are calculated with respect to strand 1. * To examine individual') 519 FORMAT(' strand 2 bases w.r.t. strand 2, reverse signs of Tip and * Tilt.') 520 FORMAT(' For Z-DNA, reverse signs of Incl and X Dsp. Y Dsp is cor *rect as printed.') 560 FORMAT(' ROLL and TILT are the simple components of base pair norm *als along minor') 561 FORMAT(' and major axes of base pairs. They are the values that *were calculated') 562 FORMAT(' in NEWHEL90 and earlier.') 563 FORMAT(' VALL, VTIL, VROL, VSLI and VTWI are the total angle betw *een base pair') 564 FORMAT(' normal vectors, and the Tilt, Roll, Slide and Twist calc *ulated relative') 565 FORMAT(' to a set of local axes halfway between each of the long *axes, the short') 566 FORMAT(' axes, and normal vectors for the two base pairs. They a *re completely') 567 FORMAT(' independent of the choice of overall helix axis.') 580 FORMAT(1X,' ANGLES BETWEEN NORMAL VECTORS TO BASE PAIRS') 581 FORMAT(1X,'I=',I3,2X,50I4) 582 FORMAT(4X,'J=',2X,50I4) C........1.........2.........3.........4.........5.........6.........712 521 FORMAT(1X,7(4X,'-',2X),2X,2F7.2,4(4X,'-',2X),4F7.2) 522 FORMAT(1X,5(3X,'_____'),2X,10(1X,'_____')) 523 FORMAT(1X,5F8.2,2X,10F6.2,' AV') 524 FORMAT(1X,5F8.2,2X,10F6.2,' SD') 555 FORMAT(//,' BASE PAIR NORMAL COSINES AND ANGLES',/) 556 FORMAT(' COS(AX) COS(AY) COS(AZ) ANG X ANG Y * ANG Z',/) C SUPPRESS LISTING OF BASE ATOM NUMBERS VIA C BEFORE TEN WRITE COMMANDS C WRITE(LUW,506) C WRITE(LUW,500) C WRITE(LUW,502)TITL C WRITE(LUW,501) NATM,NPLA REWIND L12 2020 FORMAT(3F10.5,11X,A1,I2,5A1) 2021 FORMAT(3F10.5,11X,A1,I3,4A1) DO 10 I=1,NATM IF (IFRES.LE.1) X READ (L12,2020) (ATM(J,I),J=1,3),NTP(I),IRES(I),(NAT(J,I),J=1,5) IF (IFRES.GT.1) X READ (L12,2021) (ATM(J,I),J=1,3),NTP(I),IRES(I),(NAT(J,I),J=1,4) 10 CONTINUE DO 20 I=1,50 DO 20 J=1,3 DO 20 K=1,2 20 THE(I,J,K)=0.0 IF (INPBRL.EQ.0) CALL SETBRL DO 220 M=1,NPLA C WRITE(LUW,503) M C WRITE(LUW,504) NC(1,M),NC(2,M) NA11=NB1(1,M) NA12=NB1(2,M) NA21=NB2(1,M) NA22=NB2(2,M) C WRITE(LUW,505) NA11,NA12 C WRITE(LUW,505) NA21,NA22 L=0 DO 30 IDR=NA11,NA12 L=L+1 NA1(L)=IDR 30 CONTINUE NA1(L+1)=0 C WRITE(LUW,505) (NA1(J),J=1,L) L=0 DO 40 IDR=NA21,NA22 L=L+1 NA2(L)=IDR 40 CONTINUE NA2(L+1)=0 C WRITE(LUW,505) (NA2(J),J=1,L) DO 70 K=1,15 NTEST=NA1(K) IF(NTEST) 80,80,50 50 DO 60 J=1,3 XMIN(K,J)=ATM(J,NTEST) 60 CONTINUE 70 CONTINUE K=K+1 80 NPL=K-1 CALL MATRIX(XMIN,A,NPL) CALL EIGEN2(A,EIG) DO 90 J=1,3 PLA(M,1,J)=EIG(20,J) 90 CONTINUE DO 120 K=1,15 NTEST=NA2(K) IF(NTEST) 130,130,100 100 DO 110 J=1,3 XMIN(K,J)=ATM(J,NTEST) 110 CONTINUE 120 CONTINUE K=K+1 130 NPL=K-1 CALL MATRIX(XMIN,A,NPL) CALL EIGEN2(A,EIG) DO 140 J=1,3 PLA(M,2,J)=EIG(20,J) 140 CONTINUE DO 160 K=1,15 NTEST=NA1(K) IF(NTEST.LE.0) GO TO 170 DO 150 J=1,3 XMIN(K,J)=ATM(J,NTEST) 150 CONTINUE 160 CONTINUE K=K+1 170 K=K-1 DO 190 L=1,15 NTEST=NA2(L) IF(NTEST.LE.0) GO TO 200 DO 180 J=1,3 LI=K+L XMIN(LI,J)=ATM(J,NTEST) 180 CONTINUE 190 CONTINUE L=L+1 200 NPL=K+L-1 CALL MATRIX(XMIN,A,NPL) CALL EIGEN2(A,EIG) DO 210 J=1,3 PLA(M,3,J)=EIG(20,J) 210 CONTINUE 220 CONTINUE WRITE(LUW,506) WRITE(LUW,500) WRITE(LUW,502)TITL DO 250 J=1,3 IF(J-3)600,601,600 601 WRITE(LUW,555) WRITE(LUW,556) GO TO 602 600 CONTINUE WRITE(LUW,507) J WRITE(LUW,556) 602 CONTINUE DO 240 M=1,NPLA DO 230 K=1,3 CSN=PLA(M,J,K) ANG(K)=57.29578* ACOS(CSN) 230 CONTINUE WRITE(LUW,508) (PLA(M,J,K),K=1,3),(ANG(K),K=1,3) 240 CONTINUE 250 CONTINUE C CALCULATE AND STORE TIP AND INCLINATION DO 300 I=1,NPLA NC1=NC(1,I) NC2=NC(2,I) DELX=ATM(1,NC1)-ATM(1,NC2) DELY=ATM(2,NC1)-ATM(2,NC2) DENO=SQRT(DELX**2+DELY**2) SUMX=(ATM(1,NC1)+ATM(1,NC2))*0.5 SUMY=(ATM(2,NC1)+ATM(2,NC2))*0.5 HVEC=SQRT(SUMX**2+SUMY**2) SLI(I)=(DELX*SUMX+DELY*SUMY)/DENO DIS(I)=(SUMX*DELY-DELX*SUMY)/DENO HYP(I)=HVEC ROT12(I)=DENO DELZ=ATM(3,NC1)-ATM(3,NC2) SV(1)=DELX SV(2)=DELY SV(3)=DELZ DO 260 K=1,3 RV(1,K)=PLA(I,1,K) RV(2,K)=PLA(I,2,K) 260 CONTINUE SV2=SV(1)**2+SV(2)**2+SV(3)**2 DO 280 J=1,2 RDS=RV(J,1)*SV(1)+RV(J,2)*SV(2)+RV(J,3)*SV(3) AKO=RDS/SV2 CALB=AKO*SQRT(SV2) ALB(J)=57.29578* ACOS(CALB) DO 270 K=1,3 PV(J,K)=RV(J,K)-AKO*SV(K) 270 CONTINUE PM(J)=SQRT(PV(J,1)**2+PV(J,2)**2+PV(J,3)**2) 280 CONTINUE BKL(I)=ALB(2) - ALB(1) COPT=PV(1,1)*PV(2,1)+PV(1,2)*PV(2,2)+PV(1,3)*PV(2,3) COPU=COPT/(PM(1)*PM(2)) ANSUP=57.29578* ACOS(COPU) PTW(I)=-ANSUP C REPLACE PTW(I) AS JUST CALCULATED WITH IMPROVED PROPELLER TWIST C WITH PROPER SIGN CONTROL DO 700 J=1,2 TV(J,1)=SV(2)*RV(J,3)-SV(3)*RV(J,2) TV(J,2)=SV(3)*RV(J,1)-SV(1)*RV(J,3) TV(J,3)=SV(1)*RV(J,2)-SV(2)*RV(J,1) TV2(J)=TV(J,1)**2+TV(J,2)**2+TV(J,3)**2 700 CONTINUE VSR=SV(1)*(RV(2,2)*RV(1,3)-RV(2,3)*RV(1,2))+SV(2)*(RV(2,3)*RV(1,1) *-RV(2,1)*RV(1,3))+SV(3)*(RV(2,1)*RV(1,2)-RV(2,2)*RV(1,1)) SPI=VSR*SQRT(SV2/(TV2(1)*TV2(2))) PTW(I)=57.29578*ASIN(SPI) TIL(I)=57.29578*ATAN2(DELZ,DENO) DO 290 J=1,3 SINR=(PLA(I,J,1)*DELY-PLA(I,J,2)*DELX)/DENO PHI(I,J,1)=57.29578* ASIN(SINR) SINT=(PLA(I,J,1)*DELX+PLA(I,J,2)*DELY)/DENO PHI(I,J,2)=-57.29578* ASIN(SINT) 290 CONTINUE 300 CONTINUE C CALCULATE CUP(I), AS DIFFERENCE IN SUCCESSIVE BUCKLES NST=NPLA-1 DO 610 I=1,NST I1 = I+1 CUP(I)=BKL(I1)-BKL(I) 610 CONTINUE C C C CALCULATE AND STORE ROLL AND TILT AS THETA(I,J,1) AND THETA (I,J,2) NST=NPLA-1 DO 3330 I=1,NST I1=I+1 NC1=NC(1,I) NC2=NC(2,I) ND1=NC(1,I1) ND2=NC(2,I1) DX=ATM(1,ND1)-ATM(1,ND2)-ATM(1,NC1)+ATM(1,NC2) DY=ATM(2,ND1)-ATM(2,ND2)-ATM(2,NC1)+ATM(2,NC2) DENO=SQRT(DX**2+DY**2) DO 3320 J=1,3 DA=PLA(I1,J,1)-PLA(I,J,1) DB=PLA(I1,J,2)-PLA(I,J,2) C SINR=(DA*DX+DB*DY)/DENO IF(SINR.GT.1.0) SINR=1.0 IF(SINR.LT.-1.0) SINR=-1.0 THE(I,J,1)=-57.29578*ASIN(SINR) C SINT=(DB*DX-DA*DY)/DENO IF(SINT.GT.1.0) SINT=1.0 IF(SINT.LT.-1.0) SINT=-1.0 THE(I,J,2)=57.29578*ASIN(SINT) C C C CALCULATE SLIDE DO 310 K=1,3 AM(K)=(ATM(K,NC1)-ATM(K,NC2)+ATM(K,ND1)-ATM(K,ND2))/2. BM(K)=(ATM(K,NC1)+ATM(K,NC2)-ATM(K,ND1)-ATM(K,ND2))/2. 310 CONTINUE AMAG=SQRT(AM(1)**2+AM(2)**2+AM(3)**2) ADB=AM(1)*BM(1)+AM(2)*BM(2)+AM(3)*BM(3) SLJ(I)=-ADB/AMAG SLJ(I1)=0.0 C 3320 CONTINUE 3330 CONTINUE C C C CALCULATE VECTOR FUNCTIONS VALL, VTIL, VROL, VSLI, VTWI C WRITE(LUW,506) C WRITE(LUW,1002) C1002 FORMAT(/,' CHECK NORMALITY OF BASE STEP REFERENCE VECTORS') C WRITE(LUW,516) NST=NPLA-1 DO 4000 I=1,NST I1=I+1 NC1=NC(1,I) NC2=NC(2,I) ND1=NC(1,I1) ND2=NC(2,I1) C C GET VPX, VPY, VPZ COMPONENTS OF UNIT NORMAL VECTOR VP C VPX1=PLA(I,3,1) VPY1=PLA(I,3,2) VPZ1=PLA(I,3,3) C VPX2=PLA(I1,3,1) VPY2=PLA(I1,3,2) VPZ2=PLA(I1,3,3) C C CALCULATE ZALL, TOTAL ANGLE BETWEEN SUCCESSIVE BASE PLANE NORMALS COSP=VPX1*VPX2+VPY1*VPY2+VPZ1*VPZ2 ZALL(I)=57.29578*ACOS(COSP) C C GET VLX, VLY, VLZ COMPONENTS OF UNIT LONG BASE PAIR AXIS VECTOR VL C DX1=ATM(1,NC1)-ATM(1,NC2) DY1=ATM(2,NC1)-ATM(2,NC2) DZ1=ATM(3,NC1)-ATM(3,NC2) C DX2=ATM(1,ND1)-ATM(1,ND2) DY2=ATM(2,ND1)-ATM(2,ND2) DZ2=ATM(3,ND1)-ATM(3,ND2) C D1=SQRT(DX1**2+DY1**2+DZ1**2) D2=SQRT(DX2**2+DY2**2+DZ2**2) C VLX1=DX1/D1 VLY1=DY1/D1 VLZ1=DZ1/D1 C VLX2=DX2/D2 VLY2=DY2/D2 VLZ2=DZ2/D2 C C GET VSX, VSY, VSZ COMPONENTS OF UNIT SHORT BASE PAIR AXIS VECTOR VS C VSX1=VLY1*VPZ1-VLZ1*VPY1 VSY1=VLZ1*VPX1-VLX1*VPZ1 VSZ1=VLX1*VPY1-VLY1*VPX1 C VSX2=VLY2*VPZ2-VLZ2*VPY2 VSY2=VLZ2*VPX2-VLX2*VPZ2 VSZ2=VLX2*VPY2-VLY2*VPX2 C C GET UNIT MEDIAN VECTORS VPM, VLM AND VSM C VPXME=(VPX1+VPX2)/2 VPYME=(VPY1+VPY2)/2 VPZME=(VPZ1+VPZ2)/2 C VLXME=(VLX1+VLX2)/2 VLYME=(VLY1+VLY2)/2 VLZME=(VLZ1+VLZ2)/2 C VSXME=(VSX1+VSX2)/2 VSYME=(VSY1+VSY2)/2 VSZME=(VSZ1+VSZ2)/2 C VPM=SQRT(VPXME**2+VPYME**2+VPZME**2) VLM=SQRT(VLXME**2+VLYME**2+VLZME**2) VSM=SQRT(VSXME**2+VSYME**2+VSZME**2) C VPXM=VPXME/VPM VPYM=VPYME/VPM VPZM=VPZME/VPM C VLXM=VLXME/VLM VLYM=VLYME/VLM VLZM=VLZME/VLM C VSXM=VSXME/VSM VSYM=VSYME/VSM VSZM=VSZME/VSM C C RE-ORTHOGONALIZE MEDIAN VECTORS C VSXN=VLYM*VPZM-VLZM*VPYM VSYN=VLZM*VPXM-VLXM*VPZM VSZN=VLXM*VPYM-VLYM*VPXM VSN=SQRT(VSXN**2+VSYN**2+VSZN**2) VSXM=VSXN/VSN VSYM=VSYN/VSN VSZM=VSZN/VSN C VLXN=VPYM*VSZM-VPZM*VSYM VLYN=VPZM*VSXM-VPXM*VSZM VLZN=VPXM*VSYM-VPYM*VSXM VLN=SQRT(VLXN**2+VLYN**2+VLZN**2) VLXM=VLXN/VLN VLYM=VLYN/VLN VLZM=VLZN/VLN C C C VERIFY 90 DEGREE ANGLES BETWEEN MEDIAN COORDINATE AXES C ILAB='I=' XPS='PS=' XLS='LS=' XPL='PL=' COPL=VPXM*VLXM+VPYM*VLYM+VPZM*VLZM COPS=VPXM*VSXM+VPYM*VSYM+VPZM*VSZM COLS=VSXM*VLXM+VSYM*VLYM+VSZM*VLZM ANPL=57.29578*ACOS(COPL) ANPS=57.29578*ACOS(COPS) ANLS=57.29578*ACOS(COLS) WRITE(LUW,4001)ILAB,I,XPS,ANPS,XLS,ANLS,XPL,ANPL 4001 FORMAT(1X,A2,I3,2X,A3,F8.4,2X,A3,F8.4,2X,A3,F8.4) C C C C C CALCULATE RISE, LONG AXIS SLIDE, AND SHORT AXIS SHIFT OF BASE PAIRS CX1=(ATM(1,NC1)+ATM(1,NC2))/2. CY1=(ATM(2,NC1)+ATM(2,NC2))/2. CZ1=(ATM(3,NC1)+ATM(3,NC2))/2. C CX2=(ATM(1,ND1)+ATM(1,ND2))/2. CY2=(ATM(2,ND1)+ATM(2,ND2))/2. CZ2=(ATM(3,ND1)+ATM(3,ND2))/2. C SHX=CX2-CX1 SHY=CY2-CY1 SHZ=CZ2-CZ1 C ZRIS(I)=SHX*VPXM+SHY*VPYM+SHZ*VPZM ZSLI(I)=SHX*VLXM+SHY*VLYM+SHZ*VLZM ZSLP(I)=SHX*VSXM+SHY*VSYM+SHZ*VSZM C C C Note: The simple ARCCOS will not generally give the correct sign of the C angle. The following relations is used instead: C Tan(Angle between v1 and v2) = |v1 [cross] v2| / v1 [dot] v2 C ATAN2 is the standard FORTRAN subroutine to yield the proper sign of C a tangent given as a ratio of two signed quantities. C C C CALCULATE VTWI, PROJECTED ANGLE BETWEEN BASE PAIR LONG AXES C AROUND MEAN PERPENDICULAR AXIS C C GET UNIT VECTORS VL1xVPM (OR VL1KPM), AND VL2xVPM (OR VL2KPM) VL1KPMX=VLY1*VPZM-VLZ1*VPYM VL1KPMY=VLZ1*VPXM-VLX1*VPZM VL1KPMZ=VLX1*VPYM-VLY1*VPXM C VL2KPMX=VLY2*VPZM-VLZ2*VPYM VL2KPMY=VLZ2*VPXM-VLX2*VPZM VL2KPMZ=VLX2*VPYM-VLY2*VPXM C VL1KPM=SQRT(VL1KPMX**2+VL1KPMY**2+VL1KPMZ**2) VL2KPM=SQRT(VL2KPMX**2+VL2KPMY**2+VL2KPMZ**2) C VL1KPMX=VL1KPMX/VL1KPM VL1KPMY=VL1KPMY/VL1KPM VL1KPMZ=VL1KPMZ/VL1KPM C VL2KPMX=VL2KPMX/VL2KPM VL2KPMY=VL2KPMY/VL2KPM VL2KPMZ=VL2KPMZ/VL2KPM C C FORM CROSS PRODUCT OF VL1KPM AND VL2KPM CROSSX=VL1KPMY*VL2KPMZ-VL1KPMZ*VL2KPMY CROSSY=VL1KPMZ*VL2KPMX-VL1KPMX*VL2KPMZ CROSSZ=VL1KPMX*VL2KPMY-VL1KPMY*VL2KPMX C CROSS=SQRT(CROSSX**2+CROSSY**2+CROSSZ**2) C C FORM UNIT VECTOR CROSS PRODUCT CROSSXU=CROSSX/CROSS CROSSYU=CROSSY/CROSS CROSSZU=CROSSZ/CROSS C C ESTABLISH SIGN DEPENDING ON WHETHER CROSS PRODUCT IS PARALLEL OR C ANTIPARALLEL TO THE AXIS OF ROTATION SIGNF=CROSSXU*VPXME+CROSSYU*VPYME+CROSSZU*VPZME SIMAG=SQRT(SIGNF**2) SIGN=SIGNF/SIMAG C C FIND ANGLE BETWEEN L2KPM AND L1KPM USING ATAN2 ROUTINE CROSNUM=SIGN*(SQRT(CROSSX**2+CROSSY**2+CROSSZ**2)) CROSDEN=VL2KPMX*VL1KPMX+VL2KPMY*VL1KPMY+VL2KPMZ*VL1KPMZ VTWI=57.29578*ATAN2(CROSNUM, CROSDEN) C C STORE VTWI ZTWI(I)=VTWI C C C CALCULATE VROL, PROJECTED ANGLE BETWEEN BASE PAIR NORMALS C AROUND LONG AXIS C C GET UNIT VECTORS VP1xVLM (OR P1KLM), AND VP2xVLM (OR P2KLM) P1KLMX=VPY1*VLZM-VPZ1*VLYM P1KLMY=VPZ1*VLXM-VPX1*VLZM P1KLMZ=VPX1*VLYM-VPY1*VLXM C P2KLMX=VPY2*VLZM-VPZ2*VLYM P2KLMY=VPZ2*VLXM-VPX2*VLZM P2KLMZ=VPX2*VLYM-VPY2*VLXM C P1KLM=SQRT(P1KLMX**2+P1KLMY**2+P1KLMZ**2) P2KLM=SQRT(P2KLMX**2+P2KLMY**2+P2KLMZ**2) C P1KLMX=P1KLMX/P1KLM P1KLMY=P1KLMY/P1KLM P1KLMZ=P1KLMZ/P1KLM C P2KLMX=P2KLMX/P2KLM P2KLMY=P2KLMY/P2KLM P2KLMZ=P2KLMZ/P2KLM C C C FORM CROSS PRODUCT OF P1KLM AND P2KLM CROSSX=P1KLMY*P2KLMZ-P1KLMZ*P2KLMY CROSSY=P1KLMZ*P2KLMX-P1KLMX*P2KLMZ CROSSZ=P1KLMX*P2KLMY-P1KLMY*P2KLMX C CROSS=SQRT(CROSSX**2+CROSSY**2+CROSSZ**2) C C FORM UNIT VECTOR CROSS PRODUCT OF P1KLM AND P2KLM CROSSXU=CROSSX/CROSS CROSSYU=CROSSY/CROSS CROSSZU=CROSSZ/CROSS C C C ESTABLISH SIGN DEPENDING ON WHETHER CROSS PRODUCT IS PARALLEL OR C ANTIPARALLEL TO THE AXIS OF ROTATION SIGNF=CROSSXU*VLXME+CROSSYU*VLYME+CROSSZU*VLZME SIMAG=SQRT(SIGNF**2) SIGN=SIGNF/SIMAG C C FIND ANGLE BETWEEN P2KLM AND P1KLM USING ATAN2 ROUTINE CROSNUM=SIGN*(SQRT(CROSSX**2+CROSSY**2+CROSSZ**2)) CROSDEN=P2KLMX*P1KLMX+P2KLMY*P1KLMY+P2KLMZ*P1KLMZ VROL=57.29578*ATAN2(CROSNUM, CROSDEN) C C C STORE VROL ZROL(I)=VROL C C C CALCULATE VTIL, PROJECTED ANGLE BETWEEN BASE PAIR NORMALS C AROUND SHORT AXIS C C GET UNIT VECTORS VP1xVSM (OR P1KSM), AND VP2xVSM (OR P2KSM) P1KSMX=VPY1*VSZM-VPZ1*VSYM P1KSMY=VPZ1*VSXM-VPX1*VSZM P1KSMZ=VPX1*VSYM-VPY1*VSXM C P2KSMX=VPY2*VSZM-VPZ2*VSYM P2KSMY=VPZ2*VSXM-VPX2*VSZM P2KSMZ=VPX2*VSYM-VPY2*VSXM C P1KSM=SQRT(P1KSMX**2+P1KSMY**2+P1KSMZ**2) P2KSM=SQRT(P2KSMX**2+P2KSMY**2+P2KSMZ**2) C P1KSMX=P1KSMX/P1KSM P1KSMY=P1KSMY/P1KSM P1KSMZ=P1KSMZ/P1KSM C P2KSMX=P2KSMX/P2KSM P2KSMY=P2KSMY/P2KSM P2KSMZ=P2KSMZ/P2KSM C C FORM CROSS PRODUCT OF P1KSM AND P2KSM CROSSX=P1KSMY*P2KSMZ-P1KSMZ*P2KSMY CROSSY=P1KSMZ*P2KSMX-P1KSMX*P2KSMZ CROSSZ=P1KSMX*P2KSMY-P1KSMY*P2KSMX C CROSS=SQRT(CROSSX**2+CROSSY**2+CROSSZ**2) C C FORM UNIT VECTOR CROSS PRODUCT CROSSU=CROSSX/CROSS CROSSU=CROSSY/CROSS CROSSU=CROSSZ/CROSS C C ESTABLISH SIGN DEPENDING ON WHETHER CROSS PRODUCT IS PARALLEL OR C ANTIPARALLEL TO THE AXIS OF ROTATION SIGNF=CROSSXU*VSXME+CROSSYU*VSYME+CROSSZU*VSZME SIMAG=SQRT(SIGNF**2) SIGN=SIGNF/SIMAG C C FIND ANGLE BETWEEN P2KSM AND P1KSM USING ATAN2 ROUTINE CROSNUM=SIGN*(SQRT(CROSSX**2+CROSSY**2+CROSSZ**2)) CROSDEN=P2KSMX*P1KSMX+P2KSMY*P1KSMY+P2KSMZ*P1KSMZ VTIL=57.29578*ATAN2(CROSNUM, CROSDEN) C C STORE VTIL ZTIL(I)=VTIL C 4000 CONTINUE C C CALCULATE ALL ANGLES BETWEEN NORMAL VECTORS C DO 800 I=1,NPLA DO 800 J=1,NPLA VPX1=PLA(I,3,1) VPY1=PLA(I,3,2) VPZ1=PLA(I,3,3) VPX2=PLA(J,3,1) VPY2=PLA(J,3,2) VPZ2=PLA(J,3,3) ANUM=VPX1*VPX2+VPY1*VPY2+VPZ1*VPZ2 ADE1=SQRT(VPX1**2+VPY1**2+VPZ1**2) ADE2=SQRT(VPX2**2+VPY2**2+VPZ2**2) COSP=ANUM/(ADE1*ADE2) IF(COSP.GT.1.0) COSP=1.0 IF(COSP.LT.-1.0) COSP=-1.0 ANGALL=57.29578*ACOS(COSP) IANGALL(I,J)=ANGALL 800 CONTINUE C C PRINT OUT NORMAL VECTOR ANGLE TABLE C DO 801 J=1,NPLA NJAY(J)=J 801 CONTINUE IF(NPLA.GT.25)GO TO 805 WRITE(LUW,506) WRITE(LUW,580) WRITE(LUW,502)TITL WRITE(LUW,516) WRITE(LUW,582)(NJAY(J),J=1,NPLA) WRITE(LUW,516) DO 802 I=1,NPLA WRITE(LUW,581)I,(IANGALL(I,J),J=1,NPLA) 802 CONTINUE GO TO 806 C 805 WRITE(LUW,506) WRITE(LUW,580) WRITE(LUW,502)TITL WRITE(LUW,516) WRITE(LUW,582)(NJAY(J),J=1,25) WRITE(LUW,516) DO 803 I=1,NPLA WRITE(LUW,581)I,(IANGALL(I,J),J=1,25) 803 CONTINUE WRITE(LUW,506) WRITE(LUW,580) WRITE(LUW,502)TITL WRITE(LUW,516) WRITE(LUW,582)(NJAY(J),J=26,NPLA) WRITE(LUW,516) DO 804 I=1,NPLA WRITE(LUW,581)I,(IANGALL(I,J),J=26,NPLA) 804 CONTINUE 806 CONTINUE C C PRINT OUT ROLL AND TILT ANGLES WRITE(LUW,506) WRITE(LUW,500) WRITE(LUW,502)TITL DO 335 L=1,12 SSQ(L)=0. 335 SIG(L)=0. N1=0 DO 420 J=2,3 IF(J-2)340,340,350 340 WRITE(LUW,509) WRITE(LUW,516) WRITE(LUW,510) WRITE(LUW,516) GO TO 360 350 WRITE(LUW,511) WRITE(LUW,516) WRITE(LUW,512) WRITE(LUW,516) 360 DO 410 M=1,NPLA IF(J-2)370,370,380 370 K=J-1 C........1.........2.........3.........4.........5.........6.........712 WRITE(LUW,4513)PHI(M,K,1),PHI(M,K,2),THE(M,K,1),THE(M,K,2),THE(M,K *,3),THE(M,K,4),PHI(M,J,1),PHI(M,J,2),THE(M,J,1),THE(M,J,2),THE(M,J *,3),THE(M,J,4) GO TO 400 380 N1=N1+1 SIG(1) = SIG(1) + PHI(M,J,1) SSQ(1) = SSQ(1) + PHI(M,J,1)*PHI(M,J,1) SIG(2) = SIG(2) + PHI(M,J,2) SSQ(2) = SSQ(2) + PHI(M,J,2)*PHI(M,J,2) SIG(8) = SIG(8) + PTW(M) SSQ(8) = SSQ(8) + PTW(M)*PTW(M) SIG(9) = SIG(9) + BKL(M) SSQ(9) = SSQ(9) + BKL(M)*BKL(M) SIG(11) = SIG(11) + DIS(M) SSQ(11) = SSQ(11) + DIS(M)*DIS(M) SIG(12) = SIG(12) + SLI(M) SSQ(12) = SSQ(12) + SLI(M)*SLI(M) IF (M.NE.NPLA) GO TO 390 WRITE(LUW,521)PHI(M,J,2),PHI(M,J,1),PTW(M),BKL(M),DIS(M),SLI(M) GO TO 400 390 SIG(3) = SIG(3) + THE(M,J,1) SSQ(3) = SSQ(3) + THE(M,J,1)*THE(M,J,1) SIG(4) = SIG(4) + THE(M,J,2) SSQ(4) = SSQ(4) + THE(M,J,2)*THE(M,J,2) SIG(5) = SIG(5) + THE(M,J,3) SSQ(5) = SSQ(5) + THE(M,J,3)*THE(M,J,3) SIG(6) = SIG(6) + THE(M,J,4) SSQ(6) = SSQ(6) + THE(M,J,4)*THE(M,J,4) SIG(7) = SIG(7) + TIL(M) SSQ(7) = SSQ(7) + TIL(M)*TIL(M) SIG(10) = SIG(10) + SLJ(M) SSQ(10) = SSQ(10) + SLJ(M)*SLJ(M) C........1.........2.........3.........4.........5.........6.........712 WRITE(LUW,513)ZALL(M),ZTIL(M),ZROL(M),ZSLI(M),ZTWI(M),ZRIS(M),ZSLP *(M),PHI(M,J,2),PHI(M,J,1),THE(M,J,2),THE(M,J,1),SLJ(M),CUP(M),PTW( *M),BKL(M),DIS(M),SLI(M) 400 CONTINUE 410 CONTINUE 420 CONTINUE C WRITE(LUW,522) DO 440 I=1,12 NN = N1 C IF (I.EQ.3.OR.I.EQ.4.OR.I.EQ.8) NN=N1-1 IF(I.EQ.3.OR.I.EQ.4.OR.I.EQ.5.OR.I.EQ.6.OR.I.EQ.7.OR.I.EQ.10)NN=N1 *-1 AV(I)=SIG(I)/NN SDZ = (NN*SSQ(I) - SIG(I)*SIG(I))/(NN*(NN-1)) IF (SDZ.LE.0.0001) SDZ=0. SD(I) = SQRT(SDZ) 440 CONTINUE c WRITE (LUW,523) (AV(I),I=1,12) WRITE(LUW,516) c WRITE (LUW,524) (SD(I),I=1,12) WRITE(LUW,516) WRITE(LUW,517) WRITE(LUW,518) WRITE(LUW,519) WRITE(LUW,520) WRITE(LUW,560) WRITE(LUW,561) WRITE(LUW,562) WRITE(LUW,563) WRITE(LUW,564) WRITE(LUW,565) WRITE(LUW,566) WRITE(LUW,567) RETURN END C------------------------- C SUBROUTINE CHSIGN C------------------------- SUBROUTINE CHSIGN (ICHSGN) CHARACTER * 8 ABUF CHARACTER * 4 MH,NAME CHARACTER * 1 IR,TITL CHARACTER * 3 KRES1,KRES2 CHARACTER * 3 ANAM1,ANAM2 COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),ICORL X ,L15 COMMON/ATOM/NATOM,ATM(3,2000),NAME(2,2000) COMMON/BROOK/IBARF,IFRES I=0 10 CONTINUE I=I+1 IF (I.EQ.NATOM) GO TO 100 J=I IF (IFRES.GT.1) GO TO 15 KRES1=NAME(1,I)(2:3)//' ' ANAM1=NAME(1,I)(4:4)//NAME(2,I)(1:2) GO TO 20 15 KRES1=NAME(1,I)(2:4) ANAM1=NAME(2,I)(1:3) 20 CONTINUE J=J+1 IF (J.GT.NATOM) GO TO 10 IF (IFRES.LE.1) KRES2=NAME(1,J)(2:3)//' ' IF (IFRES.GT.1) KRES2=NAME(1,J)(2:4) IF (KRES1.EQ.KRES2) GO TO 20 IF (IFRES.LE.1) ANAM2=NAME(1,J)(4:4)//NAME(2,J)(1:2) IF (IFRES.GT.1) ANAM2=NAME(2,J)(1:3) IF (ANAM1.NE.ANAM2) GO TO 20 IF (ATM(3,I).EQ.ATM(3,J)) GO TO 20 IF (ATM(3,I).GT.ATM(3,J)) RETURN ICHSGN=1 RETURN 100 WRITE (LUW,500) 500 FORMAT (//,' ***** NO MATCHING ATOMS BETWEEN RESIDUES *****',//) RETURN END C------------------------- C SUBROUTINE CYLIN C------------------------- SUBROUTINE CYLIN CHARACTER*80 INP CHARACTER * 4 MH,NAME CHARACTER * 9 NNAME(50) CHARACTER * 3 PBL,C1P,FFC,LN1,LN9,N1N9,O4P,O1P CHARACTER * 4 NTP CHARACTER * 1 IR,TITL,L1,L2,L3,NAT,CC,TT,MM,DOT,CHECK CHARACTER * 1 ZZ,XX,YY,UU,VV COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),ICORL COMMON/ATOM/NATM,ATM(3,2000),NTP(2000),NAT(5,2000) COMMON/CYLINP/INPCYL,NP(2,50),NC(2,50),NO(2,50) COMMON/BROOK/IBARF,IFRES,NPLA COMMON/INPIO/INP DIMENSION PS(2,50,5),CS(2,50,5),DQHP(2,50,4),STR(2,50,11),NK(4), 1 D1(3),D2(3),HOL(50),PDIT(50,50),AM(3),BM(3),NA(3),NB(3), 2 SIG1(7),SSQ1(7),AV1(7),SD1(7),SIG2(11),SSQ2(11),AV2(11),SD2(11) 3 ,ND(2,50),ALAMD1(50),ALAMD2(50),O4S(2,51,5),O4DIT(51,51),MP(100) 4 ,MC(100),MO(100),MD(100) C........1.........2.........3.........4.........5.........6.........712 DATA PBL/'P '/,C1P/'C1'''/,LN1/'N1 '/,LN9/'N9 '/ DATA ZZ/'Z'/,XX/'X'/,YY/'Y'/,UU/'U'/,VV/'V'/ DATA CC/'C'/,TT/'T'/,MM/'M'/,DOT/'.'/,O4P/'O4'''/,O1P/'O1'''/ 500 FORMAT(1X,' CYLINDRICAL PARAMETERS, NEWHEL93') 501 FORMAT(/,I5,' ATOMS AND ',I5,' BASE PAIRS') 502 FORMAT(/,1X,120A1,/) 503 FORMAT(2I5) 504 FORMAT(/,' PHOSPHATE BACKBONE TABLE, 5" TO 3" DIRECTION IN EACH ST *RAND',/) 505 FORMAT(' R PHI Z D Q H PI') 506 FORMAT(/,' STRAND ',I4,/) 507 FORMAT(7F8.2) 508 FORMAT(20I5) 509 FORMAT(/,2X,' ROTATION AND RISE TABLE, 5" TO 3" DIRECTION',/) 510 FORMAT(' S5" S3" Dxy TwC1" TWIST Dz RISE SLIDE *X DSP Y DSP Dxyz') 511 FORMAT(/,' STRAND ',I4,/) 512 FORMAT(5F7.2,2F7.3,4F7.2) C........1.........2.........3.........4.........5.........6.........712 513 FORMAT(' Twist, Rise, Dxyz and Dxy are magnitudes. All other quan *tities have') 514 FORMAT(' opposite signs for ascending and descending strands.') 515 FORMAT(1H1) 516 FORMAT(' PHOSPHORUS-PHOSPHORUS DISTANCES LESS 5.8 A',/) 536 FORMAT(' SUGAR O4"--O4" DISTANCES LESS 2.8 A',/) 517 FORMAT(' STRAND 1 DOWN, STRAND 2 ACROSS, 5" TO 3" DIRECTION') 555 FORMAT(1H ) 518 FORMAT(1X,20F6.2) 519 FORMAT(3F8.2,4(5X,'-',2X)) 520 FORMAT(7(2X,'______'),/) 521 FORMAT(7F8.2,' AV',/) 522 FORMAT(7F8.2,' SD',/) 523 FORMAT(F7.2,7(4X,'-',2X),3F7.2) 524 FORMAT(11(2X,'_____'),/) 525 FORMAT(11F7.2,' AV') 526 FORMAT(11F7.2,' SD',/) 527 FORMAT(4X,'-',2X,4F7.2,2F7.3,4F7.2) 528 FORMAT (//,2X,'RESIDUE LAMBDA(1) LAMBDA(2)',/) 529 FORMAT (2X,A9,2F10.2) 530 FORMAT (/,2X,'ANGLES BETWEEN GLYCOSIDIC BONDS AND C1'' - C1'' BASE X - PAIR VECTORS') RAD=0.0174533 WRITE(LUW,515) WRITE (LUW,500) WRITE(LUW,502)TITL WRITE(LUW,501) NATM,NPLA REWIND L12 2020 FORMAT(3F10.5,11X,A3,5A1) 2021 FORMAT(3F10.5,11X,A4,4A1) DO 10 I=1,NATM IF (IFRES.LE.1) X READ (L12,2020) (ATM(J,I),J=1,3),NTP(I),(NAT(J,I),J=1,5) IF (IFRES.GT.1) X READ (L12,2021) (ATM(J,I),J=1,3),NTP(I),(NAT(J,I),J=1,4) 10 CONTINUE NPHO=NPLA-1 DO 40 I=1,2 DO 30 J=1,NPHO DO 20 K=1,4 DQHP(I,J,K)=0.0 20 CONTINUE 30 CONTINUE 40 CONTINUE DO 70 I=1,2 DO 60 J=1,NPLA DO 50 K=1,11 STR(I,J,K)=0.0 50 CONTINUE 60 CONTINUE 70 CONTINUE IF (INPCYL.NE.0) GO TO 110 IFLAG=0 J=0 L=0 N=0 LO=0 DO 100 I=1,NATM N1N9=LN9 CHECK=NTP(I)(1:1) IF (CHECK.EQ.CC.OR.CHECK.EQ.TT.OR.CHECK.EQ.MM.OR.CHECK.EQ.ZZ) x N1N9=LN1 IF (CHECK.EQ.XX.OR.CHECK.EQ.YY.OR.CHECK.EQ.UU.OR.CHECK.EQ.VV) x N1N9=LN1 L1=NAT(1,I) L2=NAT(2,I) L3=NAT(3,I) FFC=L1//L2//L3 80 IF (FFC.EQ.PBL) GO TO 90 IF (FFC.EQ.O4P) GO TO 701 702 IF (FFC.NE.C1P) GO TO 95 L=L+1 MC(L)=I GO TO 100 90 J=J+1 MP(J)=I GO TO 100 701 LO=LO+1 MO(LO)=I GO TO 702 95 IF (FFC.NE.N1N9) GO TO 100 N=N+1 MD(N)=I IFLAG=1 100 CONTINUE DO 1004 M=1,NPLA NC(1,M)=MC(M) NO(1,M)=MO(M) ND(1,M)=MD(M) MR=NPLA+M NC(2,M)=MC(MR) NO(2,M)=MO(MR) ND(2,M)=MD(MR) 1004 CONTINUE DO 1005 M=1,NPHO NP(1,M)=MP(M) MR=NPHO+M NP(2,M)=MP(MR) 1005 CONTINUE C C WRITE(LUW,1001) C1001 FORMAT(/,' CHECK OF O4", C1", N1/9, P ATOM NUMBERS, STR. 1 & 2') C DO 1006 M=1,NPLA C WRITE(LUW,555) C WRITE(LUW,1007)M,NO(1,M),NO(2,M) C1007 FORMAT(I5,5X,'NO=',2I5) C WRITE(LUW,1008)M,NC(1,M),NC(2,M) C1008 FORMAT(I5,5X,'NC=',2I5) C WRITE(LUW,1009)M,ND(1,M),ND(2,M) C1009 FORMAT(I5,5X,'ND=',2I5) C WRITE(LUW,1010)M,NP(1,M),NP(2,M) C1010 FORMAT(I5,5X,'NP=',2I5) C1006 CONTINUE C WRITE(LUW,515) IF (IFLAG.EQ.0) ND(2,N)=I 110 DO 140 I=1,2 DO 130 J=1,NPHO NPH=NP(I,J) DO 120 K=1,3 PS(I,J,K)=ATM(K,NPH) 120 CONTINUE X=PS(I,J,1) Y=PS(I,J,2) R=SQRT(X**2+Y**2) PHI=57.29578*ATAN2(Y,X) PS(I,J,4)=R IF (PHI.LT.0.) PHI=360.+PHI PS(I,J,5)=PHI 130 CONTINUE 140 CONTINUE DO 141 I=1,2 DO 131 J=1,NPLA NO4=NO(I,J) DO 121 K=1,3 O4S(I,J,K)=ATM(K,NO4) 121 CONTINUE 131 CONTINUE 141 CONTINUE DO 160 I=1,2 DO 150 J=1,NPLA NCA=NC(I,J) CS(I,J,1)=ATM(1,NCA) CS(I,J,2)=ATM(2,NCA) CS(I,J,3)=ATM(3,NCA) X=CS(I,J,1) Y=CS(I,J,2) R=SQRT(X**2+Y**2) PHI=57.29578*ATAN2(Y,X) CS(I,J,4)=R CS(I,J,5)=PHI 150 CONTINUE 160 CONTINUE NPS=NPHO-1 DO 180 I=1,2 DO 170 J=1,NPS J1=J+1 DX=PS(I,J1,1)-PS(I,J,1) DY=PS(I,J1,2)-PS(I,J,2) DZ=PS(I,J1,3)-PS(I,J,3) D=SQRT(DX**2+DY**2+DZ**2) Q=SQRT(DX**2+DY**2) H=DZ SINP=H/D PIT=57.29578* ASIN(SINP) DQHP(I,J,1)=D DQHP(I,J,2)=Q DQHP(I,J,3)=H DQHP(I,J,4)=PIT 170 CONTINUE 180 CONTINUE WRITE(LUW,504) WRITE(LUW,505) DO 190 I=1,7 SSQ1(I)=0. 190 SIG1(I)=0. N1=0 DO 230 I=1,2 WRITE(LUW,506)I DO 220 J=1,NPHO N1=N1+1 SIG1(1) = SIG1(1) + PS(I,J,4) SSQ1(1) = SSQ1(1) + PS(I,J,4)*PS(I,J,4) SIG1(2) = SIG1(2) + PS(I,J,5) SSQ1(2) = SSQ1(2) + PS(I,J,5)*PS(I,J,5) SIG1(3) = SIG1(3) + PS(I,J,3) SSQ1(3) = SSQ1(3) + PS(I,J,3)*PS(I,J,3) IF (J.NE.NPHO) GO TO 200 WRITE(LUW,519) PS(I,J,4),PS(I,J,5),PS(I,J,3) GO TO 220 200 DO 210 K=1,4 SIG1(3+K) = SIG1(3+K) + DQHP(I,J,K) 210 SSQ1(3+K) = SSQ1(3+K) + DQHP(I,J,K)*DQHP(I,J,K) WRITE(LUW,507)PS(I,J,4),PS(I,J,5),PS(I,J,3),(DQHP(I,J,K),K=1,4) 220 CONTINUE 230 CONTINUE DO 240 I=1,7 NN = N1 IF (I.GT.3) NN=N1-2 AV1(I)=SIG1(I)/NN SDZ = (NN*SSQ1(I) - SIG1(I)*SIG1(I))/(NN*(NN-1)) IF (SDZ.LE.0.0001) SDZ=0. SD1(I) = SQRT(SDZ) 240 CONTINUE WRITE (LUW,520) WRITE (LUW,521) (AV1(I),I=1,7) WRITE (LUW,522) (SD1(I),I=1,7) DO 390 I=1,2 ANUL=0.000 STR(I,1,1)=ANUL STR(I,1,3)=ANUL DO 250 L=2,7 STR(I,NPLA,K)=ANUL 250 CONTINUE DO 380 J=1,NPHO J1=J+1 IJ=NPLA-J+1 IJ1=NPLA-J C CALCULATE S5' DIF=CS(I,J1,5)-PS(I,J,5) 260 IF(DIF-180.) 280,280,270 270 DIF=DIF-360. GO TO 260 280 IF(DIF+180.) 290,290,300 290 DIF=DIF+360. GO TO 280 300 STR(I,J1,1)=DIF C CALCULATE S3' DIF=PS(I,J,5)-CS(I,J,5) 310 IF(DIF-180.) 330,330,320 320 DIF=DIF-360. GO TO 310 330 IF(DIF+180.) 340,340,350 340 DIF=DIF+360. GO TO 330 350 STR(I,J,2)=DIF C CALCULATE Q(C1") AND T DX=CS(I,J,1)-CS(I,J1,1) DY=CS(I,J,2)-CS(I,J1,2) STR(I,J,3)=SQRT(DX**2+DY**2) 370 STR(I,J,4)=STR(I,J1,1)+STR(I,J,2) C CALCULATE TGLOBAL DX1=CS(1,J,1)-CS(2,IJ,1) DX2=CS(1,J1,1)-CS(2,IJ1,1) DY1=CS(1,J,2)-CS(2,IJ,2) DY2=CS(1,J1,2)-CS(2,IJ1,2) R1=SQRT(DX1**2+DY1**2) R2=SQRT(DX2**2+DY2**2) R12=DX1*DX2+DY1*DY2 COA=R12/(R1*R2) STR(I,J,5)=57.29578* ACOS(COA) C CALCULATE H AND HGLOBAL STR(I,J,6)=CS(I,J1,3)-CS(I,J,3) STR(I,J,7)=ABS(CS(1,J1,3)-CS(1,J,3)+CS(2,IJ1,3)-CS(2,IJ,3))/2. 380 CONTINUE 390 CONTINUE C PRINT ANGLE TABLE DO 420 K=5,7,2 DO 400 J=1,NPHO IJ=NPHO-J+1 HOL(J)=STR(2,IJ,K) 400 CONTINUE DO 410 J=1,NPHO STR(2,J,K)=HOL(J) 410 CONTINUE 420 CONTINUE WRITE(LUW,515) WRITE(LUW,500) WRITE(LUW,502)TITL WRITE(LUW,509) WRITE(LUW,510) DO 430 I=1,NPLA N1=NC(1,I) NK(2)=N1 NK(1)=ND(1,I) I2=NPLA-I+1 N2=NC(2,I2) NK(3)=N2 NK(4)=ND(2,I2) DO 425 K=1,3 NA(K)=NK(K) NB(K)=NK(4-K+1) 425 CONTINUE CALL ANGLE (NA,TE,ATM) ALAMD1(I)=TE/RAD CALL ANGLE (NB,TE,ATM) ALAMD2(I)=TE/RAD NNAME(I)=NTP(NK(1))//DOT//NTP(NK(3)) DX=ATM(1,N1)-ATM(1,N2) DY=ATM(2,N1)-ATM(2,N2) SX=(ATM(1,N1)+ATM(1,N2))/2. SY=(ATM(2,N1)+ATM(2,N2))/2. HM=SQRT(SX**2+SY**2) RM=SQRT(DX**2+DY**2) SL=(SX*DX+SY*DY)/RM DS=(SX*DY-SY*DX)/RM STR(1,I,10)=SL STR(1,I,9)=DS STR(1,I,11)=RM C REPLACE C1/C1 BY D(C1") DCSQ=STR(1,I,3)**2+STR(1,I,6)**2 STR(1,I,11)=SQRT(DCSQ) DCSQ=STR(2,I,3)**2+STR(2,I,6)**2 STR(2,I,11)=SQRT(DCSQ) 430 CONTINUE NST=NPLA-1 DO 450 I=1,NST I1=I+1 I2=NPLA-I+1 I3=NPLA-I NC1=NC(1,I) NC2=NC(2,I2) ND1=NC(1,I1) ND2=NC(2,I3) DO 440 K=1,3 AM(K)=(ATM(K,NC1)-ATM(K,NC2)+ATM(K,ND1)-ATM(K,ND2))/2. BM(K)=(ATM(K,NC1)+ATM(K,NC2)-ATM(K,ND1)-ATM(K,ND2))/2. 440 CONTINUE AMAG=SQRT(AM(1)**2+AM(2)**2+AM(3)**2) ADB=AM(1)*BM(1)+AM(2)*BM(2)+AM(3)*BM(3) STR(1,I,8)=-ADB/AMAG STR(1,I1,8)=0.00 450 CONTINUE DO 460 I=1,11 SSQ2(I)=0. 460 SIG2(I)=0. N1=0 DO 620 J=1,NPLA N1=N1+1 DO 470 K=9,11 SIG2(K) = SIG2(K) + STR(1,J,K) 470 SSQ2(K) = SSQ2(K) + STR(1,J,K)*STR(1,J,K) KK=2 IF (J.GT.1) GO TO 480 KK=4 SIG2(2) = SIG2(2) + STR(1,J,2) SSQ2(2) = SSQ2(2) + STR(1,J,2)*STR(1,J,2) GO TO 600 480 SIG2(1) = SIG2(1) + STR(1,J,1) SSQ2(1) = SSQ2(1) + STR(1,J,1)*STR(1,J,1) IF (J.EQ.NPLA) GO TO 620 600 DO 610 K=KK,8 SIG2(K) = SIG2(K) + STR(1,J,K) 610 SSQ2(K) = SSQ2(K) + STR(1,J,K)*STR(1,J,K) 620 CONTINUE DO 630 I=1,11 IF (I.EQ.3) THEN NN=N1-(I-1) ELSE IF (I.GE.9) THEN NN=N1 ELSE NN=N1-1 ENDIF ENDIF AV2(I)=SIG2(I)/NN SDZ = (NN*SSQ2(I) - SIG2(I)*SIG2(I))/(NN*(NN-1)) IF (SDZ.LE.0.0001) SDZ=0. SD2(I) = SQRT(SDZ) 630 CONTINUE DO 650 I=1,2 WRITE(LUW,511)I DO 640 J=1,NPLA IF (I.NE.1.OR.(J.NE.NPLA.AND.J.NE.1)) WRITE(LUW,512) * (STR(I,J,K),K=1,11) IF (I.EQ.1.AND.J.EQ.NPLA) WRITE (LUW,523) STR(I,J,1), * (STR(I,J,K),K=9,11) IF (I.EQ.1.AND.J.EQ.1) WRITE (LUW,527) STR(I,J,2), * (STR(I,J,K),K=3,11) 640 CONTINUE IF (I.GT.1) GO TO 650 WRITE(LUW,524) WRITE (LUW,525) (AV2(K),K=1,11) WRITE(LUW,508) WRITE (LUW,526) (SD2(K),K=1,11) 650 CONTINUE WRITE(LUW,555) WRITE(LUW,513) WRITE(LUW,514) WRITE(LUW,515) WRITE(LUW,500) WRITE(LUW,502)TITL WRITE(LUW,530) WRITE(LUW,528) WRITE (LUW,529) (NNAME(I),ALAMD1(I),ALAMD2(I),I=1,NPLA) C CALCULATE AND PRINT ALL P-P DISTANCES BETWEEN STRANDS WRITE(LUW,515) OPEN(UNIT=7,FILE=INP(1:(index(INP,' ')-1))//'.po4', * STATUS='UNKNOWN') WRITE(7,555) WRITE(7,555) WRITE(7,555) WRITE(7,555) WRITE(7,555) WRITE(LUW,500) WRITE(7,500) WRITE(LUW,502)TITL WRITE(7,502)TITL DO 670 I=1,NPHO DO 670 J=1,NPHO SUMMA=0.0 DO 660 K=1,3 660 SUMMA=SUMMA+(PS(1,I,K)-PS(2,J,K))**2 PDIT(I,J)=SQRT(SUMMA)-5.8 670 CONTINUE WRITE(LUW,516) WRITE(7,516) WRITE(LUW,517) WRITE(7,517) DO 680 I=1,NPHO WRITE(LUW,555) WRITE(7,555) WRITE(LUW,518)(PDIT(I,J),J=1,NPHO) WRITE(7,518)(PDIT(I,J),J=1,NPHO) 680 CONTINUE C CALCULATE AND PRINT ALL O4"-O4" DISTANCES BETWEEN STRANDS WRITE(LUW,515) WRITE(7,515) WRITE(7,555) WRITE(7,555) WRITE(7,555) WRITE(7,555) WRITE(7,555) WRITE(LUW,500) WRITE(7,500) WRITE(LUW,502) TITL WRITE(7,502) TITL NO4=NPHO+1 DO 770 I=1,NO4 DO 770 J=1,NO4 SUMMO=0.0 DO 760 K=1,3 760 SUMMO=SUMMO+(O4S(1,I,K)-O4S(2,J,K))**2 O4DIT(I,J)=SQRT(SUMMO)-2.8 770 CONTINUE WRITE(LUW,536) WRITE(7,536) WRITE(LUW,517) WRITE(7,517) DO 780 I=1,NO4 WRITE(LUW,555) WRITE(7,555) WRITE(LUW,518)(O4DIT(I,J),J=1,NO4) WRITE(7,518)(O4DIT(I,J),J=1,NO4) 780 CONTINUE RETURN END C------------------------- C SUBROUTINE DEFROT C------------------------- SUBROUTINE DEFROT(VECTOR,ANGLE,ROTOR) C DEFINE A GENERAL ROTATION MATRIX DIMENSION VECTOR(3),DIRCOS(3),ROTOR(3,3) DO 10 I=1,3 DO 10 J=1,3 10 ROTOR(I,J)=0.0 C=COS(ANGLE) S=SIN(ANGLE) VMAG=SQRT(DOTP(VECTOR,VECTOR)) DO 20 I=1,3 20 DIRCOS(I)=VECTOR(I)/VMAG ROTOR(1,2)=DIRCOS(3)*S ROTOR(2,1)=-ROTOR(1,2) ROTOR(3,1)=DIRCOS(2)*S ROTOR(1,3)=-ROTOR(3,1) ROTOR(2,3)=DIRCOS(1)*S ROTOR(3,2)=-ROTOR(2,3) DO 30 I=1,3 DO 30 J=1,3 30 ROTOR(I,J)=ROTOR(I,J)+DIRCOS(I)*DIRCOS(J)*(1.0-C) DO 40 I=1,3 40 ROTOR(I,I)=ROTOR(I,I)+C RETURN END C------------------------- C FUNCTION DOTP C------------------------- FUNCTION DOTP(X,Y) C TAKE THE DOT PRODUCT OF TWO CARTESIAN VECTORS DIMENSION X(3),Y(3) DOTP=0.0 DO 10 I=1,3 DOTP=DOTP+X(I)*Y(I) 10 CONTINUE RETURN END C------------------------- C SUBROUTINE EIGEN1 C------------------------- SUBROUTINE EIGEN1(A,R,N,MV) C COMPUTE EIGENVALUES AND EIGENVECTORS OF A REAL SYMMETRIC MATRIX DIMENSION A(1),R(1) RANGE=1.0E-6 IF(MV-1) 10,25,10 10 IQ=-N DO 20 J=1,N IQ=IQ+N DO 20 I=1,N IJ=IQ+I R(IJ)=0.0 IF(I-J) 20,15,20 15 R(IJ)=1.0 20 CONTINUE 25 ANORM=0.0 DO 35 I=1,N DO 35 J=I,N IF(I-J) 30,35,30 30 IA=I+(J*J-J)/2 ANORM=ANORM+A(IA)*A(IA) 35 CONTINUE IF(ANORM) 165,165,40 40 ANORM=1.414*SQRT(ANORM) ANRMX=ANORM*RANGE/FLOAT(N) IND=0 THR=ANORM 45 THR=THR/FLOAT(N) 50 L=1 55 M=L+1 60 MQ=(M*M-M)/2 LQ=(L*L-L)/2 LM=L+MQ 62 IF(ABS(A(LM))-THR) 130,65,65 65 IND=1 LL=L+LQ MM=M+MQ X=0.5*(A(LL)-A(MM)) 68 Y=-A(LM)/SQRT(A(LM)*A(LM)+X*X) IF(X) 70,75,75 70 Y=-Y 75 SINX=Y/SQRT(2.0*(1.0+(SQRT(1.0-Y*Y)))) SINX2=SINX*SINX 78 COSX=SQRT(1.0-SINX2) COSX2=COSX*COSX SINCS=SINX*COSX ILQ=N*(L-1) IMQ=N*(M-1) DO 125 I=1,N IQ=(I*I-I)/2 IF(I-L) 80,115,80 80 IF(I-M) 85,115,90 85 IM=I+MQ GO TO 95 90 IM=M+IQ 95 IF(I-L) 100,105,105 100 IL=I+LQ GO TO 110 105 IL=L+IQ 110 X=A(IL)*COSX-A(IM)*SINX A(IM)=A(IL)*SINX+A(IM)*COSX A(IL)=X 115 IF(MV-1) 120,125,120 120 ILR=ILQ+I IMR=IMQ+I X=R(ILR)*COSX-R(IMR)*SINX R(IMR)=R(ILR)*SINX+R(IMR)*COSX R(ILR)=X 125 CONTINUE X=2.0*A(LM)*SINCS Y=A(LL)*COSX2+A(MM)*SINX2-X X=A(LL)*SINX2+A(MM)*COSX2+X A(LM)=(A(LL)-A(MM))*SINCS+A(LM)*(COSX2-SINX2) A(LL)=Y A(MM)=X 130 IF(M-N) 135,140,135 135 M=M+1 GO TO 60 140 IF(L-(N-1)) 145,150,145 145 L=L+1 GO TO 55 150 IF(IND-1) 160,155,160 155 IND=0 GO TO 50 160 IF(THR-ANRMX) 165,165,45 165 IQ=-N DO 185 I=1,N IQ=IQ+N LL=I+(I*I-I)/2 JQ=N*(I-2) DO 185 J=I,N JQ=JQ+N MM=J+(J*J-J)/2 IF(A(LL)-A(MM)) 170,185,185 170 X=A(LL) A(LL)=A(MM) A(MM)=X IF(MV-1) 175,185,175 175 DO 180 K=1,N ILR=IQ+K IMR=JQ+K X=R(ILR) R(ILR)=R(IMR) 180 R(IMR)=X 185 CONTINUE RETURN END C------------------------- C SUBROUTINE EIGEN2 C------------------------- SUBROUTINE EIGEN2(A,EIG) CHARACTER * 1 IR,TITL COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),ICORL DIMENSION A(3,3),EIG(20,3) DIMENSION AX(6,6),B(3,3) C FIND THE LEAST EIGENVECTOR OF EQUATION AE=LE, WHERE A IS A 3X3 C MATRIX, L IS THE EIGENVALUE, AND E IS THE EIGENVECTOR. START C WITH E OR EIG=(0,0,1) AND ITERATE BY THE SWMB METHOD C (ACTA CRYST 14:185,1959) BY MULTIPLYING EIG BY B, WHERE C IS THE ADJOINT MATRIX OF A. HALT WHEN EIG CONVERGES (LESS THAN C 0.1% CHANGE BETWEEN CYCLES), OR AFTER 20 ITERATIONS. WRITE OUT C EACH EIGENVECTOR FOR INSPECTION. DO 1 I=1,3 DO 1 J=1,3 I3=I+3 J3=J+3 AX(I,J)=A(I,J) AX(I3,J)=A(I,J) AX(I,J3)=A(I,J) AX(I3,J3)=A(I,J) 1 CONTINUE DO 2 I=1,3 DO 2 J=1,3 I1=I+1 I2=I+2 J1=J+1 J2=J+2 B(J,I)=AX(I1,J1)*AX(I2,J2)-AX(I1,J2)*AX(I2,J1) 2 CONTINUE DO 3 I=2,20 DO 3 J=1,3 EIG(I,J)=0.0 3 CONTINUE EIG(1,1)=0.0 EIG(1,2)=0.0 EIG(1,3)=1.0 C SUPPRESS LISTING OF EIGENVECTOR ITERATIONS VIA C BEFORE FIVE WRITE C OR FORMAT COMMANDS C WRITE(LUW,5)(EIG(1,J),J=1,3) C5 FORMAT(' EIGENVECTOR=',3F10.5) DO 4 K=1,19 K1=K+1 EX=B(1,1)*EIG(K,1)+B(1,2)*EIG(K,2)+B(1,3)*EIG(K,3) EY=B(2,1)*EIG(K,1)+B(2,2)*EIG(K,2)+B(2,3)*EIG(K,3) EZ=B(3,1)*EIG(K,1)+B(3,2)*EIG(K,2)+B(3,3)*EIG(K,3) EN=SQRT(EX**2+EY**2+EZ**2) EIG(K1,1)=EX/EN EIG(K1,2)=EY/EN EIG(K1,3)=EZ/EN C WRITE(LUW,5)(EIG(K1,J),J=1,3) DO 6 J=1,3 DIF=ABS(EIG(K1,J)-EIG(K,J)) ARR=DIF/EIG(K1,J) ERR=ABS(ARR) IF(ERR-0.001)6,6,4 6 CONTINUE EIGZ=EIG(K1,3) SIGN=EIGZ/ABS(EIGZ) DO 8 J=1,3 EIG(20,J)=EIG(K1,J)*SIGN 8 CONTINUE GO TO 9 4 CONTINUE C TEST THE EIGENVECTOR BY CARRYING OUT AE=LE AND EXAMINING L 9 EX=A(1,1)*EIG(20,1)+A(1,2)*EIG(20,2)+A(1,3)*EIG(20,3) EY=A(2,1)*EIG(20,1)+A(2,2)*EIG(20,2)+A(2,3)*EIG(20,3) EZ=A(3,1)*EIG(20,1)+A(3,2)*EIG(20,2)+A(3,3)*EIG(20,3) ELAMX=EX/EIG(20,1) ELAMY=EY/EIG(20,2) ELAMZ=EZ/EIG(20,3) C WRITE(LUW,10)ELAMX,ELAMY,ELAMZ C10 FORMAT(' LEAST EIGENVALUES= ',3F12.8) RETURN END C------------------------- C SUBROUTINE GETATM C------------------------- SUBROUTINE GETATM(OMAT) CHARACTER * 5 ATAB CHARACTER * 4 MH,NORIG,NAME,NOX,NOY,NOZ,IBLNK CHARACTER * 1 IR,TITL,Y,R,A,T,G,C,L1,P1,P2,P3,P4,P5 CHARACTER * 60 FMT CHARACTER * 4 CC,EC COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),ICORL COMMON/INPUT/CELL(6),FMT,NPMIN,NPMAX,IFPUN,IFLAGP X ,IHELIX,IBROLL,ICYLIN,ITORNG COMMON/HELIXP/NSEQ(200,10),MMH(10),NVEC,INPHEL,NN1(200),NN2(200), X ATAB(30,10),LSEQ(10) COMMON/ATOM/NATOM,ATM(3,2000),NAME(2,2000) COMMON/BROOK/IBARF,IFRES COMMON NAT,X(3,200),VEC(3,200),PIFAC DIMENSION IW1(3),IW2(3),OMAT(3,3),OATM(3,2000),DDUM(3) DATA NORIG/'ORGN'/,NOX/'A AX'/,NOY/'B AX'/,NOZ/'C AX'/ DATA IBLNK/' '/,Y/'Y'/,R/'R'/,T/'T'/,A/'A'/,G/'G'/,C/'C'/ 1000 FORMAT(120A1) WRITE(LUW,2000) TITL 2000 FORMAT(1H1,1X,120A1,/) WRITE(LUW,2222) 2222 FORMAT(' INPUT FRAC.COORDS. X,Y,Z AND ATOM NAME IN FORMAT:') WRITE(LUW,2986)FMT 2986 FORMAT(1X,A60,/) 1010 FORMAT(2I5) IF (INPHEL.NE.0) GO TO 100 DO 90 JJ=1,IHELIX NNH=0 DO 80 I=1,NATOM L1=NAME(1,I)(1:1) C L4=NAME(1,I)(4:4) C L5=NAME(2,I)(1:1) C L6=NAME(2,I)(2:2) C L7=NAME(2,I)(3:3) IF (IFRES.LE.1) CC=NAME(1,I)(4:4)//NAME(2,I)(1:3) IF (IFRES.GT.1) CC=NAME(2,I)(1:4) C CC=L4//L5//L6//L7 NUMBR=LSEQ(JJ) DO 60 J=1,NUMBR P1=ATAB(J,JJ)(1:1) P2=ATAB(J,JJ)(2:2) P3=ATAB(J,JJ)(3:3) P4=ATAB(J,JJ)(4:4) P5=ATAB(J,JJ)(5:5) EC=P2//P3//P4//P5 IF (P1.NE.R.AND.P1.NE.Y) EC=P1//P2//P3//P4 IF (P1.EQ.R) GO TO 40 IF (P1.EQ.Y) GO TO 50 30 IF (CC.EQ.EC) GO TO 70 GO TO 60 40 IF (L1.EQ.G.OR.L1.EQ.A) GO TO 30 GO TO 60 50 IF (L1.EQ.T.OR.L1.EQ.C) GO TO 30 60 CONTINUE GO TO 80 70 NNH=NNH+1 NSEQ(NNH,JJ)=I 80 CONTINUE MMH(JJ)=NNH 90 CONTINUE 100 WRITE(LUW,2010) NATOM,IFLAGP 2010 FORMAT(' WILL READ ',I5,' ATOMS: PRINT FLAG=',I5/) NAT=NATOM ATM(1,NAT+1)=0.0 ATM(2,NAT+1)=0.0 ATM(3,NAT+1)=0.0 ATM(1,NAT+2)=1.0 ATM(2,NAT+2)=0.0 ATM(3,NAT+2)=0.0 ATM(1,NAT+3)=0.0 ATM(2,NAT+3)=1.0 ATM(3,NAT+3)=0.0 ATM(1,NAT+4)=0.0 ATM(2,NAT+4)=0.0 ATM(3,NAT+4)=1.0 NAME(1,NAT+1)=NORIG NAME(1,NAT+2)=NOX NAME(1,NAT+3)=NOY NAME(1,NAT+4)=NOZ NAME(2,NAT+1)=IBLNK NAME(2,NAT+2)=IBLNK NAME(2,NAT+3)=IBLNK NAME(2,NAT+4)=IBLNK NAT=NAT+4 1020 FORMAT(6F10.5) WRITE(LUW,2333) 2333 FORMAT(/' CELL CONSTANTS A,B,C,ALPHA,BETA,GAMMA ARE:') WRITE(LUW,2020) CELL 2020 FORMAT(1X,6(F7.3,1X)/) CALL GETMAT(OMAT,CELL) WRITE(LUW,2025) IFLAGP 2025 FORMAT(' COORDINATES IN ORTHONORMAL SYSTEM: PRINT FLAG=',I5/) DO 120 I=1,NAT DO 111 J=1,3 111 DDUM(J)=ATM(J,I) CALL ROTATE(OMAT,DDUM) DO 112 J=1,3 112 ATM(J,I)=DDUM(J) IF(IFLAGP.NE.1) GO TO 110 WRITE(LUW,2015) I,(ATM(J,I),J=1,3),(NAME(J,I),J=1,2) 110 CONTINUE 2015 FORMAT (1X,I5,5X,3(F9.4,1X),1X,2A4) DO 120 J=1,3 OATM(J,I)=ATM(J,I) 120 CONTINUE CALL MINV(OMAT,3,DETOMT,IW1,IW2) RETURN END C------------------------- C SUBROUTINE GETHVC C------------------------- SUBROUTINE GETHVC C GET VECTORS RELATING EQUIVALENT ATOMS CHARACTER * 5 ATAB CHARACTER * 4 NAME CHARACTER * 1 IR,TITL CHARACTER * 60 FMT COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),ICORL COMMON/HELIXP/NSEQ(200,10),MMH(10),NVEC,INPHEL,NN1(200),NN2(200), X ATAB(30,10),LNUMB(10) COMMON/INPUT/CELL(6),FMT,NPMIN,NPMAX,IFPUN,IFLAGP X ,IHELIX,IBROLL,ICYLIN,ITORNG COMMON/ATOM/NATOM,ATM(3,2000),NAME(2,2000) COMMON NAT,X(3,200),VEC(3,200),PIFAC DIMENSION FIELD(15) WRITE(LUW,2444) 2444 FORMAT(' **************************************'/) WRITE(LUW,2000) 2000 FORMAT(' ATOM PAIRS USED TO DETERMINE THE HELIX AXIS ARE:') IF (INPHEL.NE.0) GO TO 40 NVEC=0 DO 30 JJ=1,IHELIX NNH=MMH(JJ) KH=NNH/2 NMH=NNH-1 DO 20 K=1,NMH IF (K.EQ.KH) GO TO 20 IF (K.GT.KH) GO TO 10 NVEC=NVEC+1 IF(NVEC.GT.200) GO TO 120 NN1(NVEC)=NSEQ(K,JJ) NN2(NVEC)=NSEQ(K+1,JJ) GO TO 20 10 NVEC=NVEC+1 IF(NVEC.GT.200) GO TO 120 NN2(NVEC)=NSEQ(K,JJ) NN1(NVEC)=NSEQ(K+1,JJ) 20 CONTINUE 30 CONTINUE 40 DO 60 N=1,NVEC N2=NN2(N) N1=NN1(N) WRITE(LUW,2010) N1,(NAME(J,N1),J=1,2),N2,(NAME(J,N2),J=1,2) 2010 FORMAT(1X,2(I5,1X,2A4,8X),15A4) DO 50 I=1,3 X(I,N)=ATM(I,N1) VEC(I,N)=ATM(I,N2)-ATM(I,N1) 50 CONTINUE 60 CONTINUE WRITE(LUW,2444) IF(NVEC.GE.3) RETURN WRITE(LUW,2030) 2030 FORMAT(' LESS THAN 3 VECTORS: HELIX AXIS CANNOT BE DETERMINED') STOP 16 120 WRITE(LUW,2020) 2020 FORMAT(' ARRAYS EXCEEDED BY INPUT OF GT.200 VECTORS') STOP 16 END C------------------------- C SUBROUTINE GETMAT C------------------------- SUBROUTINE GETMAT(BMAT,CELCON) C GET THE MATRIX TO TRANSFORM INTO CARTESIAN COORDINATES REAL BMAT(3,3),CELCON(6),ANG(3) ANG(1)=CELCON(4) ANG(2)=CELCON(5) ANG(3)=CELCON(6) PIFAZE=3.1415926/180.0 DO 10 I=1,3 10 ANG(I)=ANG(I)*PIFAZE A1=CELCON(1) A2=CELCON(2) A3=CELCON(3) TEMP1=(COS(ANG(2))-COS(ANG(1))*COS(ANG(3))) TEMP2=SIN(ANG(1))*SIN(ANG(3)) OMG= ACOS(TEMP1/TEMP2) BMAT(1,1)=A1*SIN(ANG(3))*SIN(OMG) BMAT(1,2)=0.0 BMAT(1,3)=0.0 BMAT(2,1)=A1*COS(ANG(3)) BMAT(2,2)=A2 BMAT(2,3)=A3*COS(ANG(1)) BMAT(3,1)=A1*SIN(ANG(3))*COS(OMG) BMAT(3,2)=0.0 BMAT(3,3)=A3*SIN(ANG(1)) RETURN END C------------------------- C SUBROUTINE HELIX C------------------------- SUBROUTINE HELIX CHARACTER * 5 ATAB CHARACTER * 4 NAME,MH CHARACTER * 1 IR,TITL CHARACTER * 60 FMT COMMON/HELIXP/NSEQ(200,10),MMH(10),NVEC,INPHEL,NN1(200),NN2(200), X ATAB(30,10) COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),ICORL COMMON/INPUT/CELL(6),FMT,NPMIN,NPMAX,IFPUN,IFLAGP COMMON/ATOM/NATOM,ATM(3,2000),NAME(2,2000) COMMON/BROOK/IBARF,IFRES,NPLA COMMON NAT,X(3,200),VEC(3,200),PIFAC INTEGER LS1(200),LS2(200) REAL ROT(3,3),A(200,2),G(200),XYC(2),SGXYC(2),XYZC(3) REAL ROTOT(3,3),WORK(3,3),OMAT(3,3),OUHV(3) REAL X1(3),X2(3),UHV(3),XLS(200,200),ELS(200),UXYZ(3,3) DATA UXYZ/1.0, 3*0.0, 1.0, 3*0.0, 1.0/ PIFAC=3.1415926/180.0 C I/O UNIT ASSIGNMENTS:CARD READ,PRINTER,AND IFPUN.NE.0 CALL READAT CALL GETATM(OMAT) CALL GETHVC CALL MLSP(NVEC,VEC,UHV,DZ,SIGPL) TOL=1.0E-5 IF((ABS(UHV(1)).LE.TOL).AND.(ABS(UHV(2)).LE.TOL)) GO TO 200 ANG= ACOS(UHV(1)/SQRT(UHV(1)**2+UHV(2)**2)) IF(UHV(2).LE.0.0) ANG=-ANG CALL DEFROT(UXYZ(1,3),ANG,ROT) CALL MATMUL(ROTOT,ROT,UXYZ,WORK,3,3,3) DO 100 I=1,NAT CALL ROTATE(ROT,ATM(1,I)) 100 CONTINUE DO 110 I=1,NVEC CALL ROTATE(ROT,X(1,I)) CALL ROTATE(ROT,VEC(1,I)) 110 CONTINUE 200 CONTINUE ANG= ACOS(UHV(3)) CALL DEFROT(UXYZ(1,2),ANG,ROT) CALL MATMUL(ROTOT,ROT,ROTOT,WORK,3,3,3) DO 120 I=1,NAT CALL ROTATE(ROT,ATM(1,I)) 120 CONTINUE DO 130 I=1,NVEC CALL ROTATE(ROT,X(1,I)) CALL ROTATE(ROT,VEC(1,I)) 130 CONTINUE DO 210 I=1,NVEC A(I,1)=VEC(1,I)*2.0 A(I,2)=VEC(2,I)*2.0 G(I)=(VEC(1,I) + X(1,I))**2 + (VEC(2,I) + X(2,I))**2 G(I)=G(I)-X(1,I)**2 -X(2,I)**2 210 CONTINUE CALL LLSQ(NVEC,2,200,2,G,A,XYC,SGXYC,NSING,XLS,ELS,LS1,LS2) IF(NSING.EQ.0) GO TO 220 WRITE(LUW,2200) 2200 FORMAT(' >>>>>>>>>>> CENTER MATRIX SINGULAR') STOP 16 220 CONTINUE XYZC(1)=XYC(1) XYZC(2)=XYC(2) XYZC(3)=0.0 CALL MATMUL(XYZC,XYZC,ROTOT,WORK,1,3,3) WRITE(LUW,2220) 2220 FORMAT(' IN ORTHONORMAL COORDINATES, ') WRITE(LUW,2210) (UHV(JI),XYZC(JI),JI =1,3) 2210 FORMAT(' HELIX AXIS IN PARAMETRIC FORM:'/' X = ',F12.5,'*S + ', *F12.5/' Y = ',F12.5,'*S + ',F12.5/' Z = ',F12.5, '*S + ',F12.5/) CALL MATMUL(XYZC,OMAT,XYZC,WORK,3,3,1) CALL MATMUL(OUHV,OMAT,UHV,WORK,3,3,1) WRITE(LUW,2230) 2230 FORMAT(' IN ORIGINAL CRYSTAL COORDINATES,') WRITE(LUW,2210) (OUHV(JI),XYZC(JI),JI=1,3) DO 250 I=1,2 DO 230 J=1,NAT ATM(I,J)=ATM(I,J)-XYC(I) 230 CONTINUE DO 240 J=1,NVEC X(I,J)=X(I,J)-XYC(I) 240 CONTINUE 250 CONTINUE THETA=0.0 SGTHET=0.0 DO 310 I=1,NVEC DO 300 J=1,3 X1(J)=X(J,I) X2(J)=X(J,I)+VEC(J,I) 300 CONTINUE CALL POLAR(X1) CALL POLAR(X2) T21=X2(2)-X1(2) T21P=T21+PIFAC*360.0 T21M=T21-PIFAC*360.0 IF(ABS(T21P).LT.ABS(T21)) T21=T21P IF(ABS(T21M).LT.ABS(T21)) T21=T21M THETA=THETA+T21 SGTHET=SGTHET+T21*T21 310 CONTINUE THETA=THETA/FLOAT(NVEC) SGTHET=(SGTHET-FLOAT(NVEC)*THETA*THETA)/FLOAT(NVEC-1) C SET SGTHET TO ZERO IF NEGATIVE (ONLY ARISES WITH IDEAL COORDS) IF(SGTHET)450,451,451 450 WRITE(LUW,452)SGTHET 452 FORMAT(/,1X,'****WARNING: SGTHET=',F10.5,3X,'SET TO ZERO****',/) SGTHET=0.0 451 CONTINUE SGTHET=SQRT(SGTHET)/PIFAC SIG=0.0 IF(NVEC.LE.3) GO TO 340 DO 330 I=1,NVEC DO 320 J=1,3 X1(J)=X(J,I) X2(J)=X(J,I)+VEC(J,I) 320 CONTINUE CALL POLAR(X1) CALL POLAR(X2) SIG=SIG+(X1(3)+DZ-X2(3))**2 + X1(1)**2 + X2(1)**2 *-2.0*X1(1)*X2(1)*COS(X1(2)+THETA-X2(2)) 330 CONTINUE SIG=SQRT(SIG/FLOAT(NVEC-3)) 340 CONTINUE IF(DZ.LT.0.0) THETA=-THETA IF(DZ.LT.0.0) DZ=-DZ ANG=THETA/PIFAC WRITE(LUW,2300) ANG,DZ,SIG 2300 FORMAT(' >>>>>> HELIX ROTATION: ',F8.3,5X,'DISPLACEMENT: ', *F8.4/' ',10X,'STATISTICS:'/' ',20X,'OVERALL STANDARD DEV.: ', *F8.4) WRITE(LUW,2310) SGXYC,SIGPL,SGTHET 2310 FORMAT(21X,'SIGMA(X): ',F7.4,', SIGMA(Y): ',F7.4, *', SIGMA(Z)=SIGMA(DISPLACEMENT): ',F7.4,', SIGMA(ROTATION): ', *F7.3) PTURN=ABS(360.0/ANG) WRITE(LUW,2320)PTURN 2320 FORMAT(11X,'THERE ARE ',F5.2,' RESIDUES PER TURN'/) ICHSGN=0 CALL CHSIGN (ICHSGN) CALL PUTATM(THETA,DZ,ICHSGN) WRITE(LUW,2400) 2400 FORMAT(' NORMAL END OF JOB: YOU HAVE GIVEN BIRTH TO A HELIX') WRITE(LUW,2401) 2401 FORMAT(' NUMBER BASE PAIRS =') WRITE(LUW,2402)NPLA 2402 FORMAT(5X,I5) RETURN END C------------------------- C SUBROUTINE INIT C------------------------- SUBROUTINE INIT CHARACTER * 5 ATAB CHARACTER * 4 NAME,MH CHARACTER * 1 IR,TITL,IBL CHARACTER * 60 FMT COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),ICORL X ,L15 COMMON/INPUT/CELL(6),FMT,NPMIN,NPMAX,IFPUN,IFLAGP X ,IHELIX,IBROLL,ICYLIN,ITORNG COMMON/BROOK/IBARF,IFRES,NBP COMMON/ATOM/NATOM,ATM(3,2000),NAME(2,2000) COMMON/HELIXP/NSEQ(200,10),MMH(10),NVEC,INPHEL,NN1(200),NN2(200), X ATAB(30,10),LSEQ(10) COMMON/CYLINP/INPCYL,NP(2,50),NC(2,50),NO(2,50) COMMON/BROLLP/INPBRL,NB(2,50),NB1(2,50),NB2(2,50) COMMON/TORANP/INPTOR,NCHI(2,100,4),NATO(2,1000) DATA IBL/' '/ NATOM=0 ICORL=0 IFLAGP=0 IFPUN=0 NPMAX=0 NPMIN=0 INPHEL=0 INPBRL=0 INPCYL=0 INPTOR=0 IHELIX=0 IBROLL=0 ICYLIN=0 ITORNG=0 DO 10 I=1,120 10 TITL(I)=IBL NBP=0 DO 20 JJ=1,10 DO 20 I=1,30 20 ATAB(I,JJ)=' ' DO 30 I=1,2 DO 30 J=1,1000 30 NATO(I,J)=0 DO 40 I=1,2 DO 40 J=1,100 DO 40 L=1,4 40 NCHI(I,J,L)=0 RETURN END C------------------------- C SUBROUTINE LLSQ C------------------------- SUBROUTINE LLSQ(NOBS,NPARM,NOD,NPD,G,A,P,SIGP,NSING,FM,E,L1,L2) C GENERAL LINEAR LEAST-SQUARES ROUTINE REAL FM(NPARM,NPARM),E(NPARM) INTEGER L1(NPARM),L2(NPARM) DIMENSION G(NOD),A(NOD,NPD),P(NPD),SIGP(NPD) IF((NOBS.GT.NOD).OR.(NPARM.GT.NPD)) STOP 16 DO 10 I=1,NPARM E(I)=0.0 DO 10 J=1,NPARM FM(I,J)=0.0 10 CONTINUE DO 20 I=1,NPARM DO 20 J=1,NOBS E(I)=E(I)+A(J,I)*G(J) DO 20 K=1,NPARM FM(I,K)=FM(I,K)+A(J,I)*A(J,K) 20 CONTINUE CALL MINV(FM,NPARM,DET,L1,L2) NSING=0 IF(DET.NE.0.0) GO TO 30 NSING=1 RETURN 30 CONTINUE DO 40 I=1,NPARM P(I)=0.0 40 CONTINUE DO 50 I=1,NPARM DO 50 J=1,NPARM P(I)=P(I)+FM(I,J)*E(J) 50 CONTINUE SIGSQ=0.0 DO 70 I=1,NOBS DEL=G(I) DO 60 J=1,NPARM DEL=DEL-A(I,J)*P(J) 60 CONTINUE SIGSQ=SIGSQ+DEL*DEL 70 CONTINUE DO 80 I=1,NPARM SIGP(I)=SQRT(FM(I,I)*SIGSQ) 80 CONTINUE RETURN END C------------------------- C SUBROUTINE MATMUL C------------------------- SUBROUTINE MATMUL(PROD,XM,YM,WM,N,M,L) C MATRIX MULTIPLICATION DIMENSION XM(N,M),YM(M,L),PROD(N,L),WM(N,L) DO 10 I=1,N DO 10 J=1,L WM(I,J)=0.0 10 CONTINUE DO 20 I=1,N DO 20 K=1,L DO 20 J=1,M WM(I,K)=WM(I,K)+XM(I,J)*YM(J,K) 20 CONTINUE DO 30 I=1,N DO 30 J=1,L PROD(I,J)=WM(I,J) 30 CONTINUE RETURN END C------------------------- C SUBROUTINE MATRIX C------------------------- SUBROUTINE MATRIX(XMIN,A,NPL) C BUILT THE A MATRIX FROM A LIST OF ATOM COORDINATES DIMENSION XMIN(30,3),A(3,3) DIMENSION XMAJ(30,3),XAV(3) DO 1 J=1,3 XAV(J)=0.0 1 CONTINUE DO 2 K=1,NPL DO 2 J=1,3 XAV(J)=XAV(J)+XMIN(K,J) 2 CONTINUE DO 3 J=1,3 XAV(J)=XAV(J)/NPL 3 CONTINUE DO 4 K=1,NPL DO 4 J=1,3 XMAJ(K,J)=XMIN(K,J)-XAV(J) 4 CONTINUE DO 5 I=1,3 DO 5 J=1,3 A(I,J)=0.0 5 CONTINUE DO 6 K=1,NPL DO 6 I=1,3 DO 6 J=1,3 A(I,J)=A(I,J)+XMIN(K,I)*XMAJ(K,J) 6 CONTINUE RETURN END C------------------------- C SUBROUTINE MINV C------------------------- SUBROUTINE MINV(A,N,D,L,M) C INVERT A MATRIX DIMENSION A(1),L(1),M(1) D=1.0 NK=-N DO 80 K=1,N NK=NK+N L(K)=K M(K)=K KK=NK+K BIGA=A(KK) DO 20 J=K,N IZ=N*(J-1) DO 20 I=K,N IJ=IZ+I 10 IF(ABS(BIGA)-ABS(A(IJ))) 15,20,20 15 BIGA=A(IJ) L(K)=I M(K)=J 20 CONTINUE J=L(K) IF(J-K) 35,35,25 25 KI=K-N DO 30 I=1,N KI=KI+N HOLD=-A(KI) JI=KI-K+J A(KI)=A(JI) 30 A(JI)=HOLD 35 I=M(K) IF(I-K) 45,45,38 38 JP=N*(I-1) DO 40 J=1,N JK=NK+J JI=JP+J HOLD=-A(JK) A(JK)=A(JI) 40 A(JI)=HOLD 45 IF(BIGA) 48,46,48 46 D=0.0 RETURN 48 DO 55 I=1,N IF(I-K) 50,55,50 50 IK=NK+I A(IK)=A(IK)/(-BIGA) 55 CONTINUE DO 65 I=1,N IK=NK+I HOLD=A(IK) IJ=I-N DO 65 J=1,N IJ=IJ+N IF(I-K) 60,65,60 60 IF(J-K) 62,65,62 62 KJ=IJ-I+K A(IJ)=HOLD*A(KJ)+A(IJ) 65 CONTINUE KJ=K-N DO 75 J=1,N KJ=KJ+N IF(J-K) 70,75,70 70 A(KJ)=A(KJ)/BIGA 75 CONTINUE D=D*BIGA A(KK)=1.0/BIGA 80 CONTINUE K=N 100 K=(K-1) IF(K) 150,150,105 105 I=L(K) IF(I-K) 120,120,108 108 JQ=N*(K-1) JR=N*(I-1) DO 110 J=1,N JK=JQ+J HOLD=A(JK) JI=JR+J A(JK)=-A(JI) 110 A(JI)=HOLD 120 J=M(K) IF(J-K) 100,100,125 125 KI=K-N DO 130 I=1,N KI=KI+N HOLD=A(KI) JI=KI-K+J A(KI)=-A(JI) 130 A(JI)=HOLD GO TO 100 150 RETURN END C------------------------- C SUBROUTINE MLSP C------------------------- SUBROUTINE MLSP(NPOINT,X,U,D,SIG) C GENERAL LEAST-SQUARES PLANE ROUTINE DIMENSION X(3,NPOINT),U(3),XAV(3),XM(3,3),R(3,3),YM(3,3),ZM(9) EQUIVALENCE (XM(1,1),ZM(1)) IF(NPOINT.LT.3) STOP 16 DO 10 I=1,3 XAV(I)=0.0 10 CONTINUE DO 20 I=1,3 DO 20 J=1,NPOINT XAV(I)=XAV(I)+X(I,J) 20 CONTINUE DO 30 I=1,3 XAV(I)=XAV(I)/FLOAT(NPOINT) 30 CONTINUE DO 40 I=1,3 DO 40 J=1,3 XM(I,J)=0.0 YM(I,J)=0.0 40 CONTINUE DO 50 I=1,3 DO 50 J=1,3 DO 50 K=1,NPOINT YM(I,J)=YM(I,J)+(X(I,K)-XAV(I))*(X(J,K)-XAV(J)) 50 CONTINUE ZM(1)=YM(1,1) ZM(2)=YM(1,2) ZM(3)=YM(2,2) ZM(4)=YM(1,3) ZM(5)=YM(2,3) ZM(6)=YM(3,3) CALL EIGEN1(XM,R,3,0) EVAL=ZM(6) SIG=0.0 IF(NPOINT.LE.3) GO TO 60 IF(EVAL.LT.0.0) GO TO 60 SIG=SQRT(EVAL/FLOAT(NPOINT-3)) 60 CONTINUE DO 70 I=1,3 U(I)=R(I,3) 70 CONTINUE UMAG=DOTP(U,U) DO 80 I=1,3 U(I)=U(I)/UMAG 80 CONTINUE D=0.0 DO 90 I=1,3 D=D+U(I)*XAV(I) 90 CONTINUE RETURN END C------------------------- C SUBROUTINE PDBOUT C------------------------- SUBROUTINE PDBOUT (NAME1,NAME2,J,X,Y,Z) CHARACTER * 4 NAME1,NAME2 CHARACTER*80 INP COMMON/INPIO/INP OPEN(UNIT=16,FILE=INP(1:(index(INP,' ')-1))//'.PDB', * STATUS='UNKNOWN') WRITE (16,100) J,NAME1(4:4),NAME2(1:3),NAME2(4:4),NAME1(1:1),NAME1 *(2:3),X,Y,Z 100 FORMAT ('ATOM',2X,I5,2X,A1,A3,2A1,5X,A2,4X,3F8.3,' 1.0 20.0') RETURN END C------------------------- C SUBROUTINE POLAR C------------------------- SUBROUTINE POLAR(V) C CONVERT A CARTESIAN VECTOR V INTO CYLINDRICAL POLAR DIMENSION V(3) X=V(1) Y=V(2) V(1)=SQRT(X*X+Y*Y) V(2)= ACOS(X/V(1)) IF(Y.LT.0.0) V(2)=-V(2) RETURN END C------------------------- C FUNCTION PRIANG C------------------------- FUNCTION PRIANG(ANGLE) C SET ANGLE T TO BETWEEN 0.0 AND 360.0 T=ANGLE 30 IF(T.LT.360.0) GO TO 40 T=T-360.0 GO TO 30 40 IF(T.GE.0.0) GO TO 50 T=T+360.0 GO TO 40 50 CONTINUE PRIANG=T RETURN END C------------------------- C SUBROUTINE PUTATM C------------------------- SUBROUTINE PUTATM(ANG,DZ,ICHSGN) C OUTPUT ATOMIC POSITIONS FOR VARIOUS POWERS OF THE HELIX OPERATOR CHARACTER * 4 NAME,MH CHARACTER * 1 IR,TITL CHARACTER * 60 FMT CHARACTER*80 INP COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),ICORL X ,L15 COMMON/INPUT/CELL(6),FMT,NPMIN,NPMAX,IFPUN,IFLAGP COMMON/ATOM/NATOM,ATM(3,2000),NAME(2,2000) COMMON NAT,X(3,200),VEC(3,200),PIFAC COMMON/INPIO/INP DIMENSION XA(3) DO 5 I=1,NAT CALL POLAR(ATM(1,I)) 5 CONTINUE 1000 FORMAT(3I5) NPMAX=NPMAX-NPMIN+1 IF (NPMAX.LE.0) NPMAX=1 C OUTPUT ON 'IPR' ALL POWERS OF THE HELIX OPERATOR, BOTH ODD AND EVEN DO 10 I=1,NPMAX POWER=FLOAT(I+NPMIN-1) IP=I+NPMIN-1 AD=ANG*POWER ZD=DZ*POWER WRITE(LUW,2444) 2444 FORMAT(' *****************************************************'/) IF (IFLAGP.EQ.0) GO TO 30 WRITE(LUW,2000) TITL,IP 2000 FORMAT(/,1X,120A1/' NO.',5X,' R',6X,'THETA',7X,'Z',5X, * 'NAME ',4X,'X',8X,'Y',8X,'Z',8X,'POWER = ',I3/) DO 10 J=1,NAT R=ATM(1,J) T=ATM(2,J)+AD Z=ATM(3,J)+ZD XC=R*COS(T) Y=R*SIN(T) T=PRIANG(T/PIFAC) WRITE(LUW,2010) J,R,T,Z,(NAME(K,J),K=1,2),XC,Y,Z 2010 FORMAT(1X,I5,5X,3(F8.3,1X),1X,2A4,2X,3(F8.3,1X)) 10 CONTINUE C OUTPUT ON UNIT 'IPU':ONLY EVEN POWERS OF THE HELIX OPERATOR IN C DIAMOND LIST FORMAT(RECORD LENGTH=132) 30 IF(IFPUN.EQ.0) GO TO 70 NOUT=0 REWIND L12 WRITE(LUW,2666) L12 2666 FORMAT(' THE FOLLOWING DIAMOND LIST HAS BEEN OUTPUT ON UNIT:',I3) WRITE(LUW,2555)TITL 2555 FORMAT(/, 120A1,/) DO 20 I=1,NPMAX POWER=FLOAT(I+NPMIN-1) IP=I+NPMIN-1 IPMOD=MOD(IP,2) IF(IPMOD.NE.0) GO TO 20 AD=ANG*POWER ZD=DZ*POWER DO 299 J=1,NATOM R=ATM(1,J) T=ATM(2,J)+AD Z=ATM(3,J)+ZD XC=R*COS(T) Y=R*SIN(T) T=PRIANG(T/PIFAC) NOUT=NOUT+1 IF (ICHSGN.NE.0) GO TO 40 Y=-Y Z=-Z 40 IF (IFLAGP.NE.0) WRITE(LUW,2020) XC,Y,Z,(NAME(K,J),K=1,2),NOUT OPEN(UNIT=L12,FILE=INP(1:(index(INP,' ')-1))//'_L12.dat', * STATUS='UNKNOWN') OPEN(UNIT=L15,FILE=INP(1:(index(INP,' ')-1))//'_L15.dat', * STATUS='UNKNOWN') WRITE(L12,2020) XC,Y,Z,(NAME(K,J),K=1,2),NOUT WRITE (L15,2021) J,(NAME(K,J),K=1,2),XC,Y,Z CALL PDBOUT (NAME(1,J),NAME(2,J),J,XC,Y,Z) 2021 FORMAT (I5,T9,2A4,T18,'1',3F10.5,' 20.00 1.0') 2020 FORMAT (3F10.5,11X,2A4,5X,I5) 299 CONTINUE 20 CONTINUE 70 CONTINUE ENDFILE L12 ENDFILE L15 RETURN END C------------------------- C SUBROUTINE READAT C------------------------- SUBROUTINE READAT CHARACTER * 80 BUF,INP CHARACTER * 4 MH,NAME,TBRF CHARACTER * 3 CHECK,O1 CHARACTER * 1 IR,TITL,TNAME(7),IBL,FOUR,ZER,IST CHARACTER * 60 FMT COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),IREAD COMMON/INPUT/CELL(6),FMT,NPMIN,NPMAX,IFPUN,IFLAGP COMMON/ATOM/NATOM,ATM(3,2000),NAME(2,2000) COMMON/BROOK/IBARF,IFRES COMMON/INPIO/INP DATA NMAX/2000/,IBL/' '/,ZER/'0'/,O1/'O1'''/,FOUR/'4'/ DATA IST/'*'/ READ(*,'(A80)') INP OPEN(UNIT=L11,FILE=INP(1:(index(INP,' ')-1)),STATUS='OLD', * READONLY) C OPEN(UNIT=L11,STATUS='OLD',READONLY) IB=3 NOCH=6 IF (IFRES.LE.1) GO TO 5 IB=4 NOCH=7 5 MAX=NMAX IF (NATOM.NE.0) MAX=NATOM I=0 10 CONTINUE I=I+1 TNAME(7)=' ' C >>>>>>>>SET MAX. NO. OF ATOMS=ARRAY SIZE-4<<<<<<<<<<<<<<<<< IF (I.GT.MAX) GO TO 50 IF(I.LE.NMAX) GO TO 20 WRITE(LUW,500) 500 FORMAT(' THAT EXCEEDS THE ATM(I,J),NAME(I),AND OATM(I,J) ARRAYS') STOP 16 20 CONTINUE IF (IREAD.EQ.0) GO TO 120 30 CONTINUE READ(L11,501,END=50,ERR=30) BUF 501 FORMAT (A80) TBRF=BUF(1:4) IF(IBARF.EQ.2.AND.TBRF.NE.'ATOM') GO TO 30 TBRF=BUF(14:14) IF(IBARF.EQ.2.AND.TBRF.EQ.'H') GO TO 30 READ(BUF,FMT) (ATM(J,I),J=1,3),(TNAME(K),K=1,NOCH) CHECK=TNAME(IB+1)//TNAME(IB+2)//TNAME(IB+3) IF (CHECK.EQ.O1) TNAME(IB+2)=FOUR IF (TNAME(2).EQ.IBL) TNAME(2)=ZER IF (IFRES.EQ.2.AND.TNAME(3).EQ.IBL) TNAME(3)=ZER IF (TNAME(IB+3).EQ.IST) TNAME(IB+3)='''' NAME(1,I)=TNAME(1)//TNAME(2)//TNAME(3)//TNAME(4) NAME(2,I)=TNAME(5)//TNAME(6)//TNAME(7)//IBL IF(IFLAGP.NE.1) GO TO 10 WRITE(LUW,502) I,(ATM(J,I),J=1,3),(NAME(J,I),J=1,2) 502 FORMAT(1X,I5,5X,3(F9.4,1X),1X,2A4) GO TO 10 50 NATOM=I-1 RETURN 120 WRITE (LUW,505) 505 FORMAT (' **** FORMAT FOR READING ATOMS NOT SPECIFIED, USE *KONN* 1, *CORL*, *BRKH* OR *FFMT*',/) STOP END C------------------------- C SUBROUTINE ROTATE C------------------------- SUBROUTINE ROTATE(ROTOR,VECTOR) C ROTATE A VECTOR VIA THE MATRIX ROTOR DIMENSION ROTOR(3,3),VECTOR(3),WORK(3) DO 10 I=1,3 WORK(I)=0.0 10 CONTINUE DO 20 I=1,3 DO 20 J=1,3 WORK(I)=WORK(I)+ROTOR(I,J)*VECTOR(J) 20 CONTINUE DO 30 I=1,3 VECTOR(I)=WORK(I) 30 CONTINUE RETURN END C------------------------- C SUBROUTINE SETBRL C------------------------- SUBROUTINE SETBRL CHARACTER * 4 MH CHARACTER * 1 IR,TITL,NTP,NAT,A,G,T,C,P,D,II,MM,ZZ CHARACTER * 1 XX,YY,UU,VV INTEGER * 2 IRES CHARACTER * 3 C6,C8,FF,FC,FD COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),ICORL COMMON/ATOM/NATM,ATM(3,2000),NTP(2000),IRES(2000),NAT(5,2000) COMMON/BROOK/IBARF,IFRES,NPLA COMMON/BROLLP/INPBRL,NC(2,50),NB1(2,50),NB2(2,50) DIMENSION IGR(2,200),IGN(200),IDUM(200) DATA A/'A'/,G/'G'/,T/'T'/,C/'C'/,P/'P'/,D/''''/ DATA II/'I'/,MM/'M'/,ZZ/'Z'/,UU/'U'/,VV/'V'/,XX/'X'/,YY/'Y'/ DATA C6/'C6 '/,C8/'C8 '/ NOCH=5 IF (IFRES.GT.1) NOCH=4 IRN=1 IGR(1,1)=1 IGN(1)=IRES(1) DO 10 I=2,NATM IF (IRES(I-1).EQ.IRES(I)) GO TO 10 IGR(2,IRN)=I-1 IRN=IRN+1 IGR(1,IRN)=I IGN(IRN)=IRES(I) 10 CONTINUE IGR(2,IRN)=NATM IRNH=IRN/2 DO 100 L=1,IRNH IS=IGR(1,L) IF=IGR(2,L) FF=C8 IF (NTP(IS).EQ.A.OR.NTP(IS).EQ.G.OR.NTP(IS).EQ.II) GO TO 20 C IF (NTP(IS).EQ.A.OR.NTP(IS).EQ.G) GO TO 20 FF=C6 20 DO 30 I=IS,IF FC=NAT(1,I)//NAT(2,I)//NAT(3,I) IF (FC.NE.FF) GO TO 30 NC(1,L)=I 30 CONTINUE M=0 DO 50 I=IS,IF DO 40 J=1,NOCH IF (NAT(J,I).EQ.D.OR.NAT(J,I).EQ.P) GO TO 50 40 CONTINUE M=M+1 IDUM(M)=I 50 CONTINUE NB1(1,L)=IDUM(1) NB1(2,L)=IDUM(M) IS=IGR(1,IRN-L+1) IF=IGR(2,IRN-L+1) FD=C8 IF (NTP(IS).EQ.A.OR.NTP(IS).EQ.G.OR.NTP(IS).EQ.II) GO TO 200 C IF (NTP(IS).EQ.A.OR.NTP(IS).EQ.G) GO TO 200 FD=C6 200 DO 60 I=IS,IF FC=NAT(1,I)//NAT(2,I)//NAT(3,I) IF (FC.NE.FD) GO TO 60 NC(2,L)=I 60 CONTINUE M=0 DO 80 I=IS,IF DO 70 J=1,NOCH IF (NAT(J,I).EQ.D.OR.NAT(J,I).EQ.P) GO TO 80 70 CONTINUE M=M+1 IDUM(M)=I 80 CONTINUE NB2(1,L)=IDUM(1) NB2(2,L)=IDUM(M) 100 CONTINUE RETURN END C------------------------- C SUBROUTINE SETTOR C------------------------- SUBROUTINE SETTOR CHARACTER * 4 MH CHARACTER * 1 IR,TITL,NTP,NAT,A,G INTEGER * 2 IRES CHARACTER * 3 FF(6),FC(4),FD(4),ATAB1(6),ATAB2(4),ATAB3(4),CHECK X ,ATAB4(5),FR(5) COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),ICORL COMMON/ATOM/NATM,ATM(3,2000),NTP(2000),IRES(2000),NAT(5,2000) COMMON/TORANP/INPTOR,NCHI(2,100,4),NATO(2,1000),NRNG(2,100,5) COMMON/BROOK/IBARF,IFRES,NPLA DIMENSION IGR(2,200),IGN(200),IDUM(200) DATA ATAB1/'P ','O5''','C5''','C4''','C3''','O3'''/ DATA ATAB4/'C4''','O4''','C1''','C2''','C3'''/ DATA ATAB2/'C4 ','N9 ','C1''','O4'''/ DATA ATAB3/'C2 ','N1 ','C1''','O4'''/ DATA A/'A'/,G/'G'/ IRN=1 IGR(1,1)=1 IGN(1)=IRES(1) DO 10 I=2,NATM IF (IRES(I-1).EQ.IRES(I)) GO TO 10 IGR(2,IRN)=I-1 IRN=IRN+1 IGR(1,IRN)=I IGN(IRN)=IRES(I) 10 CONTINUE IGR(2,IRN)=NATM IRNH=IRN/2 C GENERATE ATOM NUMBERS FOR 1ST AND 2ND STRANDS DO 200 M=1,2 NUP=0 DO 110 K=1,IRNH L=K+(M-1)*IRNH N=0 NNR=0 IS=IGR(1,L) IF=IGR(2,L) DO 20 J=1,4 20 FD(J)=ATAB2(J) IF (NTP(IS).EQ.A.OR.NTP(IS).EQ.G) GO TO 40 DO 30 J=1,4 30 FD(J)=ATAB3(J) 40 DO 60 J=1,6 DO 50 I=IS,IF CHECK=NAT(1,I)//NAT(2,I)//NAT(3,I) IF (CHECK.NE.ATAB1(J)) GO TO 50 NUP=NUP+1 IF (NUP.GT.1000) GO TO 160 NATO(M,NUP)=I GO TO 60 50 CONTINUE 60 CONTINUE DO 80 J=1,4 DO 70 I=IS,IF CHECK=NAT(1,I)//NAT(2,I)//NAT(3,I) IF (CHECK.NE.FD(J)) GO TO 70 N=N+1 IF (N.GT.4) GO TO 150 NCHI(M,K,N)=I GO TO 80 70 CONTINUE 80 CONTINUE DO 100 J=1,5 DO 90 I=IS,IF CHECK=NAT(1,I)//NAT(2,I)//NAT(3,I) IF (CHECK.NE.ATAB4(J)) GO TO 90 NNR=NNR+1 IF (NNR.GT.5) GO TO 170 NRNG(M,K,NNR)=I GO TO 90 90 CONTINUE 100 CONTINUE 110 CONTINUE 200 CONTINUE RETURN 150 WRITE (LUW,500) FD,L 500 FORMAT (//,' ***** INPUT ERROR, NUMBER OF ATOMS OF TYPE ', X 4(2X,A3,','),' IN RESIDUE ',I4,' EXCEEDS 4') STOP 160 WRITE (LUW,501) ATAB1 501 FORMAT (//,' ***** INPUT ERROR, NUMBER OF ATOMS OF TYPE ', X 6(2X,A3,','),' EXCEEDS 1000') STOP 170 WRITE (LUW,502) ATAB4,L 502 FORMAT (//,' ***** INPUT ERROR, NUMBER OF ATOMS OF TYPE ', X 5(2X,A3,','),' IN RESIDUE ',I4,' EXCEEDS 5') STOP END C------------------------- C SUBROUTINE TORANG C------------------------- SUBROUTINE TORANG CHARACTER * 4 MH CHARACTER * 8 NNAME(2,100,5) CHARACTER * 1 IR,TITL,NTP,NATT,SLASH CHARACTER * 2 BTM(3) INTEGER * 2 IRES COMMON/BROOK/IBARF,IFRES,NBP COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),ICORL COMMON/ATOM/NAT,ATM(3,2000),NTP(2000),IRES(2000),NATT(5,2000) COMMON/TORANP/INPTOR,NCHI(2,100,4),NATO(2,1000),NRNG(2,100,5) DIMENSION NA(4),ANG(300),ACHI(100),ANSW(200,14),ARG(100,5) X ,SIG(14),SSQ(14),ANGL(2,200,5),AV(14),SD(14) Y ,AVV(2,100),THAU(2,100),FAN(200),ANSX(200) DATA SLASH/'-'/ 500 FORMAT(9F8.1) 501 FORMAT(/,' STRAND',I5) 502 FORMAT(/,' ALPHA BETA GAMMA DELTA EPSILON ZETA CHI 1 DIF EP-ZE',/) 503 FORMAT(/,' ALPHA BETA GAMMA DELTA EPSILON ZETA CHI 1 MEAN EP-ZE',/) 504 FORMAT(/,' V0 V1 V2 V3 V4 Pseud. Delta 1',/) 505 FORMAT(I5,' ATOMS AND',I5,' BASE PAIRS') 506 FORMAT(1X,' TORSION ANGLE CALCULATION, NEWHEL93') 507 FORMAT (/,1X,120A1,/) 508 FORMAT(8(2X,'______'),/) 509 FORMAT(7F8.2,8X,' AV',/) 510 FORMAT(7F8.2,8X,' SD',///) 511 FORMAT(2(5X,'-',2X),7F8.1) 512 FORMAT (4F8.1,2(5X,'-',2X),3F8.1) 513 FORMAT(6(2X,'______'),/) 514 FORMAT(6F8.2,' AV',/) 515 FORMAT(6F8.2,' SD',///) 516 FORMAT (/,' SUGAR ANGLES',//,2X,5(A8,2X),' AVERAGE THAU(M)' * ,/) 517 FORMAT (6F10.1,F12.4) 518 FORMAT(6F8.1,'(',F5.1,')') 519 FORMAT (1X,'SUGAR PSEUDOROTATION PARAMS, NEWHEL93') 555 FORMAT(1H1) C READ STARTING PARAMETERS RAD=0.0174533 DANG=SIN(36*RAD)+SIN(72*RAD) WRITE(LUW,555) WRITE(LUW,506) WRITE(LUW,507)TITL WRITE(LUW,505)NAT,NBP C READ ATOM COORDINATES 2020 FORMAT(3F10.5,11X,A1,I2,5A1) 2021 FORMAT(3F10.5,11X,A1,I3,4A1) REWIND L12 DO 10 I=1,NAT IF (IFRES.LE.1) X READ (L12,2020) (ATM(J,I),J=1,3),NTP(I),IRES(I),(NATT(J,I),J=1,5) IF (IFRES.GT.1) X READ (L12,2021) (ATM(J,I),J=1,3),NTP(I),IRES(I),(NATT(J,I),J=1,4) 10 CONTINUE C GET MAIN CHAIN AND CHI ATOM NUMBERS IF (INPTOR.EQ.0) CALL SETTOR 50 CONTINUE DO 60 JC=1,14 SIG(JC)=0. SSQ(JC)=0. DO 60 IC=1,50 60 ANSW(IC,JC)=0. N1=0 DO 65 M=1,2 DO 65 I=1,100 65 AVV(M,I)=0. DO 200 M=1,2 C CALCULATE MAIN CHAIN TORSION ANGLES ANG(1)=0.0 ANG(2)=0.0 NMA=6*NBP-4 DO 80 I=1,NMA DO 70 J=1,4 IU=I+J-1 70 NA(J)=NATO(M,IU) CALL TORSON(NA,TE,ATM) I2=I+2 80 ANG(I2)=TE NMB=NMA+3 NMC=NMA+4 ANG(NMB)=0.0 ANG(NMC)=0.0 DO 120 I=1,NBP DO 110 K=1,5 DO 90 J=1,4 IU=K+J-1 IF (IU.GT.5) IU=IU-IU/5*5 NN=NRNG(M,I,IU) 90 NA(J)=NRNG(M,I,IU) CALL TORSON(NA,TE,ATM) IF(TE-180.)669,669,670 670 TE=TE-360. 669 CONTINUE ARG(I,K)=TE DO 100 J=1,3 IU=K+J-1 C NNAME(M,I,K)=BTM(1)//SLASH//BTM(2)//SLASH//BTM(3) IF (IU.GT.5) IU=IU-IU/5*5 NN=NRNG(M,I,IU) BTM(J)=NATT(1,NN)//NATT(2,NN) 100 NA(J)=NN CALL ANGLE (NA,TE,ATM) ANGL(M,I,K) = TE/RAD NNAME(M,I,K)=BTM(1)//SLASH//BTM(2)//SLASH//BTM(3) AVV(M,I)=AVV(M,I)+ANGL(M,I,K) 110 CONTINUE 120 CONTINUE C CALCULATE GLYCOSYL CHI ANGLES DO 140 I=1,NBP DO 130 J=1,4 130 NA(J)=NCHI(M,I,J) CALL TORSON(NA,TE,ATM) 140 ACHI(I)=TE C REARRANGE ANGLES AND PRINT OUT DO 180 I=1,NBP N1=N1+1 IM=NBP*(M-1)+I DO 150 J=1,6 NO=6*I-6+J 150 ANSW(IM,J)=ANG(NO) ANSW(IM,7)=ACHI(I) DO 160 J=1,5 160 ANSW(IM,8+J)=ARG(I,J) ANSW(IM,14)=(ATAN(((ARG(I,5)+ARG(I,2))-(ARG(I,4)+ARG(I,1)))/ X (2.*ARG(I,3)*DANG)))/RAD IF(ARG(I,3))671,672,672 671 ANSW(IM,14)=ANSW(IM,14)+180. 672 CONTINUE FAN(IM)=ANSW(IM,14) IF (FAN(IM).LT.0.) FAN(IM)=360.+ANSW(IM,14) 180 CONTINUE 200 CONTINUE M=1 WRITE(LUW,501)M WRITE(LUW,502) DO 210 I=1,NBP IM=NBP*2-I+1 DEL=ANSW(I,4)-ANSW(IM,4) ANSW(I,8)=DEL ANSX(I)=ANSW(I,5)-ANSW(I,6) IF(I.GT.1.AND.I.NE.NBP) WRITE(LUW,500)(ANSW(I,J),J=1,8),ANSX(I) IF(I.EQ.1) WRITE(LUW,511)(ANSW(I,J),J=3,8),ANSX(I) IF(I.EQ.NBP) WRITE(LUW,512)(ANSW(I,J),J=1,4),(ANSW(I,J),J=7,8),ANS *X(I) 210 CONTINUE M=2 WRITE(LUW,501)M WRITE(LUW,503) DO 220 I=1,NBP IM=I+NBP IN=NBP-I+1 FMEAN=0.5*(ANSW(IM,4)+ANSW(IN,4)) ANSW(IM,8)=FMEAN ANSX(IM)=ANSW(IM,5)-ANSW(IM,6) IF(I.GT.1.AND.I.NE.NBP) WRITE(LUW,500)(ANSW(IM,J),J=1,8),ANSX(IM) IF(I.EQ.1) WRITE(LUW,511)(ANSW(IM,J),J=3,8),ANSX(IM) IF(I.EQ.NBP) WRITE(LUW,512)(ANSW(IM,J),J=1,4),(ANSW(IM,J),J=7,8),A *NSX(IM) 220 CONTINUE DO 230 M=1,2 DO 230 I=1,NBP IM=NBP*(M-1)+I DO 230 J=1,14 SIG(J) = SIG(J) + ANSW(IM,J) 230 SSQ(J) = SSQ(J) + ANSW(IM,J)*ANSW(IM,J) DO 240 I=1,14 NN = N1 IF (I.LT.3.OR.(I.GT.4.AND.I.LT.7)) NN=N1-2 AV(I)=SIG(I)/NN SDZ = (NN*SSQ(I) - SIG(I)*SIG(I))/(NN*(NN-1)) IF (SDZ.LE.0.0001) SDZ=0. SD(I) = SQRT(SDZ) 240 CONTINUE WRITE (LUW,508) WRITE (LUW,509) (AV(J),J=1,7) WRITE (LUW,510) (SD(J),J=1,7) WRITE(LUW,555) WRITE (LUW,519) WRITE(LUW,507)TITL DO 260 M=1,2 WRITE(LUW,501)M WRITE (LUW,504) DO 250 K=1,NBP I = (M-1)*NBP + K WRITE (LUW,500) (ANSW(I,J),J=9,14),ANSW(I,4) 250 CONTINUE 260 CONTINUE WRITE (LUW,513) WRITE (LUW,514) (AV(J),J=9,14) WRITE (LUW,515) (SD(J),J=9,14) WRITE(LUW,555) WRITE(LUW,507)TITL DO 280 M=1,2 WRITE(LUW,501)M WRITE (LUW,516) (NNAME(M,1,J),J=1,5) DO 270 K=1,NBP AVV(M,K)=AVV(M,K)/5 THAU(M,K)=SQRT((108.-AVV(M,K))/0.00186) WRITE (LUW,517) (ANGL(M,K,J),J=1,5),AVV(M,K),THAU(M,K) 270 CONTINUE 280 CONTINUE RETURN END C------------------------- C SUBROUTINE TORSON C------------------------- C TORSION ANGLE CALCULATING SUBROUTINE SUBROUTINE TORSON(NA,TE,ATM) DIMENSION ATM(3,2000),NA(4) DIMENSION AV(3),BV(3),CV(3),S(3),T(3) DO 1 J=1,3 J1=J+1 NS=NA(J) NE=NA(J1) AV(J)=ATM(1,NE)-ATM(1,NS) BV(J)=ATM(2,NE)-ATM(2,NS) 1 CV(J)=ATM(3,NE)-ATM(3,NS) S(1)=BV(1)*CV(2)-BV(2)*CV(1) S(2)=CV(1)*AV(2)-CV(2)*AV(1) S(3)=AV(1)*BV(2)-AV(2)*BV(1) T(1)=BV(2)*CV(3)-BV(3)*CV(2) T(2)=CV(2)*AV(3)-CV(3)*AV(2) T(3)=AV(2)*BV(3)-AV(3)*BV(2) RM=SQRT(AV(2)**2+BV(2)**2+CV(2)**2) SM=SQRT(S(1)**2+S(2)**2+S(3)**2) TM=SQRT(T(1)**2+T(2)**2+T(3)**2) TOP=AV(1)*(BV(2)*CV(3)-BV(3)*CV(2))+ *AV(2)*(BV(3)*CV(1)-BV(1)*CV(3))+ *AV(3)*(BV(1)*CV(2)-BV(2)*CV(1)) SAL=(RM*TOP)/(SM*TM) CAL=(S(1)*T(1)+S(2)*T(2)+S(3)*T(3))/(SM*TM) TE=57.29578*ATAN2(SAL,CAL) IF(TE)667,668,668 667 TE=TE+360. 668 CONTINUE RETURN END C------------------------- C SUBROUTINE XDRED C------------------------- SUBROUTINE XDRED CHARACTER * 5 ATAB CHARACTER* 4 MH,I,M,JZ,NAME CHARACTER*1 TITL,IR,IBLNK,M1 CHARACTER * 60 FMT,FMKONN(2),FMCORL,FMBRKH CHARACTER * 1 AIR(4,30,10) COMMON/ATOM/NATM,ATM(3,2000),NAME(2,2000) COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),IREAD DIMENSION AA(200) COMMON/HELIXP/NSEQ(200,10),MMH(10),NVEC,INPHEL,NN1(200),NN2(200), X ATAB(30,10),KSEQ(10) COMMON/INPUT/CELL(6),FMT,NPMIN,NPMAX,IFPUN,IFLAGP X ,IHELIX,IBROLL,ICYLIN,ITORNG COMMON/CYLINP/INPCYL,NP(2,50),NC(2,50),NO(2,25) COMMON/BROLLP/INPBRL,NBB(2,50),NB1(2,50),NB2(2,50) COMMON/TORANP/INPTOR,NCHI(2,100,4),NATO(2,1000),NRNG(2,1000) COMMON/BROOK/IBARF,IFRES,NBP DATA FMBRKH/'(T31,3F8.3,T20,A1,T25,2A1,T14,3A1)'/ DATA FMKONN/'(T19,3F10.5,T9,6A1)', X '(T19,3F10.5,T8,7A1)'/ DATA FMCORL/'(T16,3F10.5,T4,A1,T10,5A1)'/ DATA JZ/' '/,IBLNK/' '/ NA=0 IFRES=1 IREAD=0 KSEQ(1)=0 10 FORMAT(A4,120A1) 20 FORMAT(1H ,A4,120A1) C C READ IN FREE FORMAT IBARF=0 30 READ(LUR,10,END=820)I,IR NA=0 NB=0 IF(I.NE.JZ)GO TO 35 WRITE (LUW,20) I,IR GO TO 30 35 N=0 NK=37 40 IF(I.EQ.MH(NK)) GOTO 50 NK=NK-1 IF(18-NK)40,50,50 50 NK=NK-17 GO TO 180 60 W=1. 70 V=0. NB=0 Y=1. U=10. Z=1. GOTO 100 80 Z=Y*Z V=U*ABS(V)+Z*XX NB=1 IF(V)90,100,90 90 V=SIGN(V,W) W=V 100 N=N+1 K=10 IF(120-N)240,110,110 110 XX=0. DO 120 MM=1,10 M1=MH(MM)(1:1) IF(IR(N).EQ.MH(MM)) GOTO 80 120 XX=XX+1. GOTO 140 130 M1=MH(K+9)(1:1) IF(IR(N).EQ.M1) GOTO 240 140 K=K-1 IF(K-2)240,130,130 150 U=1. Y=0.1 GOTO 100 160 READ(LUR,10,END=820)M,IR WRITE(LUW,20)M,IR N=0 IF(JZ.EQ.M)GOTO 60 170 WRITE (LUW,5000) 5000 FORMAT (/,' **************** INPUT ERROR') GOTO 830 180 CONTINUE 220 DO 230 J=1,200 230 AA(J)=0. NA=0 WRITE(LUW,20)I,IR GOTO 60 240 IF(K-2) 250,150,250 250 NA=NA+NB IF(200-NA)170,260,260 260 IF(-NA)270,280,280 270 AA(NA)=V+AA(NA) 280 IF(K-9)290,160,350 290 IF(K-3)60,300,60 300 W=-1. GOTO 70 350 GOTO(490,400,410,420,430,440,450,460,470,480,500,510,520,530,540 X,550,560,570,580,820) , NK C UNITCELL 400 CONTINUE IF (NA.LT.6) GO TO 170 DO 402 J=1,6 402 CELL(J)=AA(J) GO TO 30 C FPUN 410 CONTINUE IFPUN=1 GO TO 30 C PMIN 420 CONTINUE NPMIN=AA(1) GO TO 30 C PMAX 430 CONTINUE NPMAX=AA(1) GO TO 30 C KONN 440 CONTINUE FMT = FMKONN(IFRES) IREAD=1 GO TO 30 C CORL 450 CONTINUE IREAD=2 FMT=FMCORL GO TO 30 C FFMT 460 CONTINUE DO 461 J=1,72 461 FMT(J:J)=IR(J) IREAD=2 GO TO 30 C HELX 470 CONTINUE IHELIX=IHELIX+1 NNUMBR=0 J=0 471 CONTINUE J=J+1 IF (J.GT.120) GO TO 474 IF (IR(J).EQ.IBLNK) GOTO 471 NNUMBR=NNUMBR+1 AIR(1,NNUMBR,IHELIX)=IR(J) AIR(2,NNUMBR,IHELIX)=IR(J+1) AIR(3,NNUMBR,IHELIX)=IR(J+2) AIR(4,NNUMBR,IHELIX)=IR(J+3) J=J+3 GO TO 471 474 KSEQ(IHELIX)=NNUMBR GO TO 30 C NATM 480 CONTINUE IF (AA(1).NE.0) NATM=AA(1) GO TO 30 C DUMM 490 CONTINUE GO TO 30 C TITL 500 DO 501 J=1,120 501 TITL(J)=IR(J) GO TO 30 C FLGP 510 CONTINUE IFLAGP=1 GO TO 30 C BASE 520 CONTINUE NBP=AA(1) GO TO 30 C BROL 530 CONTINUE IBROLL=1 IF (NA.EQ.0) GO TO 30 INPBRL=1 IF (NBP*6.EQ.NA) GO TO 531 WRITE (LUW,5004) 5004 FORMAT (/,' ****** INPUT ERROR, NUMBER OF INPUT PARAMETERS DOES NO XT MATCH, CHECK NUMBER OF BASE PAIRS') STOP 531 JJ=0 DO 532 J=1,NBP JJ=JJ+1 NBB(1,J)=AA(JJ) NBB(2,J)=AA(JJ+1) NB1(1,J)=AA(JJ+2) NB1(2,J)=AA(JJ+3) NB2(1,J)=AA(JJ+4) NB2(2,J)=AA(JJ+5) JJ=JJ+5 532 CONTINUE GO TO 30 C CYLN 540 CONTINUE ICYLIN=1 IF (NA.EQ.0) GO TO 30 INPCYL=1 NBPM=NBP-1 IF ((2*NBPM+2*NBP).EQ.NA) GO TO 541 WRITE (LUW,5004) STOP 541 DO 543 J=1,NBP JJ=NBPM+J IF (J.EQ.NBP) GO TO 542 NP(1,J)=AA(J) NP(2,J)=AA(JJ) 542 JJ=JJ+NBPM NC(1,J)=AA(JJ) JJ=JJ+NBP NC(2,J)=AA(JJ) 543 CONTINUE GO TO 30 C TRNG 550 CONTINUE ITORNG=1 IF (NA.EQ.0) GO TO 30 INPTOR=1 KK=NBP*6 LL=NBP*4 IF ((KK+LL).EQ.NA) GO TO 551 WRITE (LUW,5004) STOP 551 DO 552 J=1,KK 552 NATO(1,J)=AA(J) DO 553 J=1,NBP DO 553 JJ=1,4 553 NCHI(1,J,JJ)=AA(KK+(J-1)*4+JJ) GO TO 30 C HLX2 560 CONTINUE IHELIX=1 IF (NA.EQ.0) GO TO 30 INPHEL=1 NVEC=NA/2 DO 561 J=1,NVEC JJ=(J-1)*2 NN1(J)=AA(JJ+1) NN2(J)=AA(JJ+2) 561 CONTINUE GO TO 30 C BRKH 570 CONTINUE IBARF=2 IREAD=2 FMT=FMBRKH GO TO 30 C FRES 580 CONTINUE IFRES=AA(1) GO TO 30 C END 820 CONTINUE DO 828 JJ=1,IHELIX IF (KSEQ(1).LE.0) GO TO 827 NUMBR=KSEQ(JJ) DO 826 J=1,NUMBR ATAB(J,JJ)=AIR(1,J,JJ)//AIR(2,J,JJ)//AIR(3,J,JJ)//AIR(4,J,JJ) 826 CONTINUE GO TO 828 827 KSEQ(JJ)=1 ATAB(1,JJ)='C1'' ' 828 CONTINUE IF (IREAD.NE.0) RETURN IREAD=1 FMT = FMKONN(IFRES) RETURN 830 STOP END C------------------------- C BLOCK DATA C------------------------- BLOCK DATA CHARACTER* 4 MH CHARACTER* 1 IR,TITL COMMON/DATA/IR(120),MH(37),LUR,LUW,L11,L12,L13,L14,TITL(120),ICORL X ,L15 DATA MH(1)/'0'/,MH(2)/'1'/,MH(3)/'2'/,MH(4)/'3'/,MH(5)/'4'/ DATA MH(6)/'5'/,MH(7)/'6'/,MH(8)/'7'/,MH(9)/'8'/,MH(10)/'9'/ DATA MH(11)/'.'/,MH(12)/'-'/,MH(13)/'+'/,MH(14)/'X'/,MH(15)/'Y'/ DATA MH(16)/'Z'/,MH(17)/','/,MH(18)/'='/,MH(19)/'CELL'/ DATA MH(20)/'FPUN'/,MH(21)/'PMIN'/,MH(22)/'PMAX'/,MH(23)/'KONN'/ DATA MH(24)/'CORL'/,MH(25)/'FFMT'/,MH(26)/'HELX'/,MH(27)/'NATM'/ DATA MH(28)/'TITL'/,MH(29)/'FLGP'/,MH(30)/'BASE'/,MH(31)/'BROL'/ DATA MH(32)/'CYLN'/,MH(33)/'TRNG'/,MH(34)/'HLX2'/,MH(35)/'BRKH'/ DATA MH(36)/'FRES'/,MH(37)/'END '/ DATA LUR/5/,LUW/6/,L11/11/,L12/12/,L13/13/,L14/14/,L15/15/ END