CCL Home Page
Up Directory CCL puck.f
      PROGRAM PUCK
C**** JACS,97,1354(1975)
      DIMENSION X(200),Y(200),Z(200),XC(200),YC(200),ZC(200),XR(15),YR(1
     &5),ZR(15),EN(3),EM(3),EL(3),XP(15),YP(15),ZP(15),QC(10),QS(10),Q(1
     &0),PHI(10),TLE(10),NAME(200),IAT(20,15),AF(200),IN(2),NRA(20)
1     FORMAT(4X,10A4)
2     FORMAT(1X,' BAD DATA'/' UNSUITABLE ',A4,' RECORD')
3     FORMAT(1X,' INPUT DATA'//' TITL ',10A4)
4     FORMAT(I3,3X,A4,3F10.4)
5     FORMAT(1X,'PUCKERING ANALYSIS OF ',10A4)
6     FORMAT(1X,'ATOM COORDINATES')
7     FORMAT(1X,A4,5X,3F13.4,5X,3F13.4)
8     FORMAT(//,I5,' MEMBERED RING')
9     FORMAT(1X,'CARTESIAN AND PUCKERING COORDINATES')
10    FORMAT(1X,'CELL CONSTANTS',6F10.4)
11    FORMAT(1X,'PUCKERING PARAMETERS')
12    FORMAT(1X,'M =',I3,'  Q(M) =',F7.4,'  PHI(M) =',F8.4)
13    FORMAT(1X,'Q(N/2) =',F7.4)
14    FORMAT(1X,'PUCKERING AMPLITUDE =',F7.4,'  THETA =',F8.4)
15    FORMAT(1X,'ENANTIOMORPHIC COORDINATES')
      OPEN(UNIT=1,FILE='PUCIN',STATUS='UNKNOWN')
      OPEN(UNIT=3,FILE='PUCPRT',STATUS='UNKNOWN')
      READ(1,1)TLE
      WRITE(3,5)TLE
      WRITE(3,3)TLE
      NORS=0
      NAT=0
      ICART = 0
      K=1
20    CALL FREFOR(AF,IN)
      INT=IN(2)+2
      GO TO(145,26,25,21,22,23),INT
21    IF(AF(1).LT.0) GO TO 16
      IF(AF(1).EQ.1) K=0
      ICART=AF(1)
      GO TO 20
16    K=0
      GO TO 20
22    A=AF(1)
      B=AF(2)
      C=AF(3)
      GAL=AF(4)
      GBE=AF(5)
      GGA=AF(6)
      GO TO 20
23    NORS=NORS+1
      M=AF(1)
      NRA(NORS)=AF(1)
      DO 24 J=1,M
24    IAT(NORS,J)=AF(J+1)
      GO TO 20
25    NAT=NAT+1
      NAME(NAT)=IN(1)
      X(NAT)=AF(K+1)
      Y(NAT)=AF(K+2)
      Z(NAT)=AF(K+3)
      GO TO 20
26    FACT=0.01745329
      IF(ICART.EQ.2)GO TO 40
      WRITE(3,6)
      DO 27 I=1,NAT
27    WRITE(3,4)I,NAME(I),X(I),Y(I),Z(I)
      IF(ICART.NE.1)GO TO 40
      DO 30 I=1,NAT
      XC(I)=X(I)
      YC(I)=Y(I)
30    ZC(I)=Z(I)
      GO TO 60
40    IF(ICART.EQ.0)GO TO 45
      WRITE(3,15)
      XF=0.0
      YF=0.0
      ZF=0.0
      DO 41 I=1,NAT
      X(I)=-X(I)
      Y(I)=-Y(I)
      Z(I)=-Z(I)
      XF=XF+X(I)
      YF=YF+Y(I)
41    ZF=ZF+Z(I)
      NX=1-XF/NAT
      NY=1-YF/NAT
      NZ=1-ZF/NAT
      DO 43 I=1,NAT
      X(I)=X(I)+NX
      Y(I)=Y(I)+NY
      Z(I)=Z(I)+NZ
43    WRITE(3,4)I,NAME(I),X(I),Y(I),Z(I)
45    AL=GAL*FACT
      BE=GBE*FACT
      GA=GGA*FACT
      PA=COS(BE)*COS(GA)
      R=(COS(AL)-PA)/SIN(GA)
      QA=SQRT(1.-COS(BE)**2-R**2)
      CAB=B*COS(GA)
      CAC=C*COS(BE)
      CBB=B*SIN(GA)
      CBC=C*R
      CCC=C*QA
      DO 50 I=1,NAT
      XC(I)=X(I)*A+Y(I)*CAB+Z(I)*CAC
      YC(I)=Y(I)*CBB+Z(I)*CBC
50    ZC(I)=Z(I)*CCC
60    DO 140 I=1,NORS
      NRAT=NRA(I)
      WRITE(3,8)NRAT
65    DO 70 J=1,NRAT
      IF(NRAT.EQ.NAT)IAT(I,J)=J
      XR(J)=XC(IAT(I,J))
      YR(J)=YC(IAT(I,J))
70    ZR(J)=ZC(IAT(I,J))
      BX=0.0
      BY=0.0
      BZ=0.0
      DO 80 J=1,NRAT
      BX=BX+XR(J)
      BY=BY+YR(J)
80    BZ=BZ+ZR(J)
      BX=BX/NRAT
      BY=BY/NRAT
      BZ=BZ/NRAT
      XD=0.0
      YD=0.0
      ZD=0.0
      XDD=0.0
      YDD=0.0
      ZDD=0.0
      DO 90 J=1,NRAT
      XR(J)=XR(J)-BX
      YR(J)=YR(J)-BY
      ZR(J)=ZR(J)-BZ
      SINF=SIN(6.283185*(J-1)/NRAT)
      COSF=COS(6.283185*(J-1)/NRAT)
      XD=XD+XR(J)*SINF
      YD=YD+YR(J)*SINF
      ZD=ZD+ZR(J)*SINF
      XDD=XDD+XR(J)*COSF
      YDD=YDD+YR(J)*COSF
90    ZDD=ZDD+ZR(J)*COSF
      A=YD*ZDD-ZD*YDD
      B=ZD*XDD-XD*ZDD
      C=XD*YDD-YD*XDD
      R=SQRT(A**2+B**2+C**2)
      EN(1)=A/R
      EN(2)=B/R
      EN(3)=C/R
      D=EN(1)*XR(1)+EN(2)*YR(1)+EN(3)*ZR(1)
      XO=XR(1)-EN(1)*D
      YO=YR(1)-EN(2)*D
      ZO=ZR(1)-EN(3)*D
      R=SQRT(XO**2+YO**2+ZO**2)
      EM(1)=XO/R
      EM(2)=YO/R
      EM(3)=ZO/R
      EL(1) =EM(2)*EN(3)-EM(3)*EN(2)
      EL(2)=EM(3)*EN(1)-EM(1)*EN(3)
      EL(3)=EM(1)*EN(2)-EM(2)*EN(1)
      A=0.0
      B=0.0
      IE=0
      WRITE(3,9)
      DO 100 J=1,NRAT
      XP(J)=EL(1)*XR(J)+EL(2)*YR(J)+EL(3)*ZR(J)
      YP(J)=EM(1)*XR(J)+EM(2)*YR(J)+EM(3)*ZR(J)
      ZP(J)=EN(1)*XR(J)+EN(2)*YR(J)+EN(3)*ZR(J)
100   WRITE(3,7)NAME(IAT(I,J)),XC(IAT(I,J)),YC(IAT(I,J)),ZC(IAT(I,J)),XP
     &(J),YP(J),ZP(J)
      M=(NRAT-1)/2
      WRITE(3,11)
      DO 120 N=2,M
      QS(N)=0.0
      QC(N)=0.0
      DO 110 J=1,NRAT
      QC(N)=QC(N)+ZP(J)*COS(6.283185*N*(J-1)/NRAT)
110   QS(N)=QS(N)+ZP(J)*SIN(6.283185*N*(J-1)/NRAT)
      IF(NRAT/2.GT.M)IE=1
      QC(N)=QC(N)*SQRT(2./NRAT)
      QS(N)=-QS(N)*SQRT(2./NRAT)
      Q(N)=SQRT(QC(N)**2+QS(N)**2)
      SNTH=QS(N)/Q(N)
      CSTH=QC(N)/Q(N)
      IF(CSTH.LT.0)GO TO 115
      PHI(N)=ASIN(SNTH)/FACT
      IF(PHI(N).LT.0)PHI(N)=PHI(N)+360.0
      GO TO 120
115   PHI(N)=ACOS(CSTH)/FACT
      IF(SNTH.LT.0)PHI(N)=360.0-PHI(N)
120   WRITE(3,12)N,Q(N),PHI(N)
      IF(IE.NE.1)GO TO 140
      BQ=0.0
      N2=NRAT/2
      Q(N2)=0.0
      DO 130 J=1,NRAT
      Q(N2)=(-1)**(J-1)*ZP(J)+Q(N2)
130   BQ=BQ+ZP(J)**2
      Q(N2)=Q(N2)*SQRT(1./NRAT)
      WRITE(3,13)Q(N2)
      BQ=SQRT(BQ)
      TH=ACOS(Q(N2)/BQ)/FACT
      WRITE(3,14)BQ,TH
140   CONTINUE
      GO TO 150
145   WRITE(3,2)IN(1)
      CLOSE(1)
      CLOSE(3)
150   STOP
      END
      SUBROUTINE FREFOR(A,IN)
      DIMENSION IR(76),IH(14),A(200),IN(2),JZ(8)
      DATA IH(1)/1H0/,IH(2)/1H1/,IH(3)/1H2/,IH(4)/1H3/,IH(5)/1H4/
      DATA IH(6)/1H5/IH(7)/1H6/,IH(8)/1H7/,IH(9)/1H8/,IH(10)/1H9/
      DATA IH(11)/1H./,IH(12)/1H-/,IH(13)/1H+/,IH(14)/1H=/
      DATA JZ(1)/4H    /,JZ(2)/4HFVAR/,JZ(3)/4HWGHT/,JZ(4)/4HAFIX/
      DATA JZ(5)/4HCOOR/,JZ(6)/4HCELL/,JZ(7)/4HRING/,JZ(8)/4HEND /
1     FORMAT(A4,76A1)
2     FORMAT(1H ,A4,76A1)
3     READ(1,1)I,IR
      IN(1)=I
      WRITE(3,2)I,IR
      IN(2)=1
      DO 5 J=1,4
      IF(I.EQ.JZ(J))GO TO 3
5     CONTINUE
      IF(I.EQ.JZ(8)) GO TO 40
      DO 15 J=5,7
      IF(I.EQ.JZ(J))IN(2)=J-3
15    CONTINUE
      N=0
      GO TO 25
6     W=1.
7     V=0.
      NB=0
      Y=1.
      U=10.
      Z=1.
      GO TO 10
8     Z=Y*Z
      V=U*ABS(V)+Z*X
      NB=1
      IF(V)9,10,9
9     V=SIGN(V,W)
      W=V
10    N=N+1
      K=6
      IF(76-N)27,11,11
11    X=0.
      DO 12 M=1,10
      IF(IR(N).EQ.IH(M))GO TO 8
12    X=X+1.
      GO TO 14
13    IF(IR(N).EQ.IH(K+9)) GO TO 27
14    K=K-1
      IF(K-2)27,13,13
17    U=1.
      Y=0.1
      GO TO 10
21    READ(1,1)M,IR
      WRITE(3,2)M,IR
      N=0
      IF(JZ(1).EQ.M)GO TO 6
C   ERROR MESSAGE RETURN.  IN(2)=-1
23    IN(2)=-1
      GO TO 35
25    DO 26 J=1,200
26    A(J)=0.
      NA=0
      GO TO 6
27    IF(K-2)28,17,28
28    NA=NA+NB
      IF(200-NA)23,29,29
29    IF(-NA)30,31,31
30    A(NA)=V+A(NA)
31    IF(K-5)33,21,35
32    CONTINUE
33    IF(K-3)6,34,6
34    W=-1.
      GO TO 7
40    IN(2)=0
C   END OF FILE
35    RETURN
      END
Modified: Mon Apr 15 16:00:00 1996 GMT
Page accessed 7937 times since Sat Apr 17 21:34:30 1999 GMT