      SUBROUTINE COEF(J,K)
      include 'param.inc'
      COMMON/COEFF/AC(10,16),POSI(30),SCAT(30),DFP(16),DFPP(16)
      COMMON/CNV/N,SINLAM(30),F(30)
      DIMENSION AB(10)
      EXTERNAL CV1,CV2
      DO 10 I=1,K
      SINLAM(I)=POSI(I)*POSI(I)
   10 F(I)=SCAT(I)
      N=K
      NN=MIN0(N,9)
      DO 15 I=1,NN
   15 AB(I)=I
      NN1=NN+1
      DO 16 I=NN1,10
   16 AB(I)=0
      CALL STEEP(AB,12,10,CV2,CV1)
      IF(AB(3).EQ.0.)AB(3)=-1E-6
      DO 20 I=1,9
   20 AC(I,J)=AB(I)
      RETURN
      END
      SUBROUTINE CV1(AB,V)
      include 'param.inc'
      COMMON/CNV/N,SINLAM(30),F(30)
      DIMENSION AB(10),V(10)
      MAXX=MIN0(9,N)
      DO 30 K=1,MAXX,2
      V(K)=0.0
      V(K+1)=0.0
      DO 30 I=1,N
      T=F(I)
      DO 20 J=1,MAXX,2
   20 T=T-AB(J)*EXP(-AB(J+1)*SINLAM(I))
      T1=-2.*T*EXP(-AB(K+1)*SINLAM(I))
      T2=-T1*SINLAM(I)*AB(K)
      V(K)=V(K)+T1
   30 V(K+1)=V(K+1)+T2
      NN=MAXX+1
      DO 10 I=NN,10
   10 V(I)=0.0
      RETURN
      END
      SUBROUTINE CV2(AB,G)
      include 'param.inc'
      COMMON/CNV/N,SINLAM(30),F(30)
      DIMENSION AB(10)
      MAXX=MIN0(9,N)
      G=0
      DO 30 I=1,N
      T=F(I)
      DO 20 J=1,MAXX,2
   20 T=T-AB(J)*EXP(-AB(J+1)*SINLAM(I))
   30 G=G+T*T
      RETURN
      END
      SUBROUTINE REFGEN(IPHASE,ZERO)
      include 'param.inc'
      COMMON/HKLCTL/KZ(3,24),AZ(24),NCN(70),ICQ,N1,IER,HN
      COMMON/CNTRLS/JOBTYP,NPROF,NPHASE,IOT,IPL,IPC,MAT,NXT,LST1,LST2,
     *LST3,IPL1,IPL2,MCYCLE,EPS,MAXS,IFOUR,NFMT
      COMMON/JNK/ALOW(30),AHIGH(30),POS(30),BCK(30),
     *NATOM(NUPH),PREF(NUPH,3),NBCKGD,NEXCRG,NSCAT,IWT
      COMMON/REFLS/IREFS(IRS),REFS(IRS,8),MLT(IRS),ICX(NUPH)
      REAL LAMDA
      INTEGER H1,H2,H3
      COMMON/G1/THMIN,STEP,THMAX,U,V,W,LAMDA(2),TMV,CTHM,RLIM,SB,WDT,
     *UU,VV,WW,SYNPOL
      COMMON/SPG/NSPGRP,K,M(10),FILL(3)
      COMMON/CELLX/A,B,C,ALPHA,BETA,GAMMA,AL(3,3)
      EQUIVALENCE (AL(1,1),ASTAR),(AL(2,2),BSTAR),(AL(3,3),CSTAR),
     *(AL(1,2),BACOSC),(AL(1,3),CACOSB),(AL(2,3),BCCOSA)
   24 FORMAT (///21H TOO MANY REFLECTIONS)
   36 FORMAT(15H LAUE SYMMETRY ,A8,33H WILL BE USED TO GENERATE INDICES)
      CHARACTER*8 LAU(14)
      DATA LAU(1),LAU(2),LAU(3),LAU(4),LAU(5),LAU(6),LAU(7),LAU(8),
     *LAU(9),LAU(10),LAU(11),LAU(12),LAU(13),LAU(14)
     */'1BAR','2/M','MMM','4/M','4/MMM','3BAR   R','3BAR M R','3BAR',
     *'3BAR M 1','3BAR 1 M','6/M','6/MMM','M3','M3M'/
      SIN1(A)=SIN(A*6.28318531)
      RAD=3.14159265359/360.0
      IF ( K.GT.3 ) STOP 5001
      KXIS=K
      K=1
 1002 N = 0
      J = 0
      I2D = 0
      I3D = 0
      I1D = 1
      I12D = 1
      I13D = 1
      I123D = 1
      I23D = 1
      I2DEL = 1
      I3DEL = 1
      DO 8000 I=1,8
      IF ( M(I).EQ.0 ) GO TO 9000
      IF ( M(I)-8192 ) 8002,8001,8000
 8001 I3DEL = 3
      I1D = 2
      I123D = 3
      GO TO 9000
 8002 IF ( M(I)-576 ) 8004,8003,8004
 8003 I3DEL = 2
      I23D = 2
      GO TO 8000
 8004 IF ( M(I)-516 ) 8008,8009,8008
 8009 I3DEL = 2
      I13D = 2
      GO TO 8000
 8008 IF ( M(I)-68 ) 8012,8013,8012
 8012 IF ( M(I)-580 ) 8015,8016,8015
 8013 I2DEL = 2
      I12D = 2
      GO TO 8000
 8016 I123D = 2
      I3DEL = 2
 8015 CONTINUE
 8000 CONTINUE
 9000 N=J-1
 97   TMIN = (SIN1(THMIN/720.0)/LAMDA(1))**2
CY    CODE FOR REFLECTIONS BEYOND OBS PATTERN
      IF(IPHASE.GE.2) GO TO 1
      TTMAX=TAN(THMAX*RAD)
      HWMAX=SQRT(U*TTMAX*TTMAX+V*TTMAX+W)
C     HKL = 0 TO BYPASS ANISOTROPIC SIZE CALCULATION
      IF(NPROF.EQ.3.OR.NPROF.EQ.5) HWMAX=VHALF(IPHASE,TTMAX,0,0,0)
      NEXCRG=NEXCRG+1
      ALOW(NEXCRG)=THMAX + STEP
      THMAX=THMAX+HWMAX*WDT/2.0
      AHIGH(NEXCRG)=THMAX
      GO TO 2
    1 THMAXO=ALOW(NEXCRG)-STEP
      TTMAX2=TAN(THMAXO*RAD)
      HWMAX2=SQRT(U*TTMAX2*TTMAX2 + V*TTMAX2 + W)
      IF(NPROF.EQ.3.OR.NPROF.EQ.5) HWMAX2=VHALF(IPHASE,TTMAX2,0,0,0)
      THMAX2=THMAXO + HWMAX2*WDT/2.0
      IF(THMAX2.LE.THMAX) GO TO 2
      THMAX=THMAX2
      AHIGH(NEXCRG)=THMAX
    2 CONTINUE
CY    END CODE
      IC = 0
      IF(IPHASE.GE.2) THEN
      DO 3 NP=1,IPHASE-1
    3 IC=IC+ICX(NP)
      ENDIF
      SMAX = (SIN1(THMAX/720.0)/LAMDA(1))**2
      IF(LAMDA(2).GT.LAMDA(1))TMIN=TMIN*(LAMDA(1)/LAMDA(2))**2
      IF(LAMDA(2).GT.0..AND.LAMDA(2).LT.LAMDA(1))SMAX=SMAX*
     *(LAMDA(1)/LAMDA(2))**2
      CALL OP1
      WRITE(16,36)LAU(NSPGRP)
      I1MAX=A*2.0*SQRT(SMAX)
      I2MAX=B*2.0*SQRT(SMAX)
      I3MAX=C*2.0*SQRT(SMAX)
      L1 = NSPGRP
      IF ( L1.GE.13.AND.M(1).EQ.580 ) L1=L1+2
      IF ( I12D+I13D+I23D.NE.5 ) GO TO 2131
      I12D = 2
      I13D = 1
      I23D = 2
      I2DEL = 2
      I3DEL = 2
 2131 IF ( L1.EQ.6 ) I2D=1
      IF ( L1.EQ.7 ) I2D=2
      I3D = I2D
      I1 = -1
 2200 I1 = 1+I1
      IF ( I1-I1MAX ) 2201,2201,2303
 2201 H1 = I1
      IF ( I2D.GT.0 ) I2MAX=I1
      GO TO (2203,2204,2204,2205,2206,2208,2209,2207,2204,2206,2204,
     1 2206,2206,2206,2206,2206),L1
 2203 I2= -MIN0(I2MAX,I1*I2MAX)
      GO TO 2210
 2204 I2 = 0
      GO TO 2210
 2205 I2= MIN0(I1,1)
      GO TO 2210
 2206 I2 = I1
      GO TO 2210
 2207 I2= MIN0(-I1+1,0)
      GO TO 2210
 2208 I2 = -2*I1
      GO TO 2210
 2209 I2 = -I1/2
 2210 I2= I2DEL*(I2/I2DEL)+MOD(I1,I12D)
      GO TO 2221
 2220 I2 = I2DEL+I2
 2221 H2 = I2
      IF ( I2-I2MAX ) 2222,2222,2200
 2222 GO TO (2223,2223,2224,2224,2224,2225,2225,2224,2224,2223,2224,
     1 2224,2226,2227,2236,2235),L1
 2223 I3= -MIN0(I3MAX,(I1+IABS(I2))*I3MAX)
      GO TO 2228
 2224 I3 = 0
      GO TO 2228
 2225 I3 = -I1-I2
      GO TO 2228
 2226 I3= MIN0(I2,I1+I3DEL)
      GO TO 2228
 2227 I3 = I2
 2228 I3= I3DEL*(I3/I3DEL)+MOD(MOD(I1+I1D*I2,I123D)+I123D,I123D) +MOD(I1
     1,I13D)+MOD(I2,I23D)
      IF ( I3D-1 ) 2232,2229,2231
 2229 I3MAX = I1
      IF ( I2-I1 ) 2230,2232,2230
 2230 I3MAX = I3MAX-1
      GO TO 2232
 2231 I3MAX = I2
      GO TO 2232
 2235 I3= I2+MOD(I1,2)
      GO TO 2232
 2236 I3 = I1+2-MOD(I2,2)
      IF ( I1.EQ.I2.AND.MOD(I1,2).EQ.0 ) I3=I1
      GO TO 2232
 2233 I3 = I3DEL+I3
 2232 H3 = I3
      IF ( KXIS.NE.3 ) GO TO 2240
      H3 = I2
      H2 = I3
 2240 CONTINUE
      IF ( I3-I3MAX ) 113,113,2220
  113 SQ=FLOAT(H1*H1)*ASTAR+FLOAT(H2*H2)*BSTAR+FLOAT(H3*H3)*CSTAR
     1+FLOAT(H1*H2)*BACOSC+FLOAT(H1*H3)*CACOSB+FLOAT(H2*H3)*BCCOSA
      SQ=SQ/4.
 2234 IF ( SQ-SMAX ) 3000,3000,2233
 2303 ICX(IPHASE)=IC
      IF(IPHASE.GE.2) THEN
      DO 4 NP=1,IPHASE-1
    4 ICX(IPHASE)=ICX(IPHASE)-ICX(NP)
      ENDIF
      CALL SORT(IPHASE)
      RETURN
 3000 IF ( SQ - TMIN ) 2233,3117,3117
 3117 IF ( IC.GT.IRS-2 ) GO TO 6001
      KZ(1,1)=H1
      KZ(2,1)=H2
      KZ(3,1)=H3
      AZ(1)=0.0
      IER=0
      CALL SMTRY2
      IF(IER.NE.0)GOTO 2233
      LXN=2
      IF(LAMDA(2).EQ.0..OR.LAMDA(2).EQ.LAMDA(1)) LXN=1
      DO 3118 LX1=1,LXN
      SQH=SQ*LAMDA(LX1)*LAMDA(LX1)
      TAN2=SQH/(1.-SQH)
      IF(SQH.GE.1.)GOTO 3118
      TANX=SQRT(TAN2)
      PIT=ATAN(TANX)/RAD+ZERO
      IF(PIT.GT.THMAX.OR.PIT.LT.THMIN) GOTO 3118
      IC=IC+1
      IF(TMV.LE.0.000001) GOTO 3127
      ROUSE=-(1.7133-0.0368*SQH)*TMV+(0.0927+0.375*SQH)*TMV*TMV
      TRANS=EXP(ROUSE)
      GO TO 3128
 3127 TRANS=1.0
 3128 CONTINUE
C-----ALL DATA
      PLOR=1./(2.*SQH*SQRT(1.-SQH))*TRANS
C-----CHARACTERISTIC X-RAY DATA
      IF(JOBTYP.EQ.3) PLOR=PLOR*(1.+(1.-2.*SQH)*(1.-2.*SQH)*CTHM)
      IF(JOBTYP.EQ.1) PLOR=PLOR*(1.+(1.-2.*SQH)*(1.-2.*SQH)*CTHM)
C-----SYNCHROTRON DATA
      IF (SYNPOL.NE.0.0) THEN
      PLOR = 1.0-SYNPOL+SYNPOL*(1.-2.*SQH)*(1.-2.*SQH)*CTHM
      PLOR = PLOR / (SQH*SQRT(1.-SQH))*TRANS
      ENDIF
      IREFS(IC)=256*(256*(256*(4*IPHASE+LX1)+128+H1)+128+H2)+128+H3
      MLT(IC)=MULT(H1,H2,H3)
      REFS(IC,1)=SQRT(U*TAN2+V*TANX+W)
      IF(NPROF.EQ.3.OR.NPROF.EQ.5)
     +       REFS(IC,1)=VHALF(IPHASE,TANX,H1,H2,H3)
      REFS(IC,2)=PIT
      REFS(IC,3)=PLOR
      IF(NPROF.NE.3.AND.NPROF.NE.5)
     * REFS(IC,4)=WW*REFS(IC,2)*REFS(IC,2) + VV*REFS(IC,2) + UU
 3118 CONTINUE
      GOTO 2233
 6001 WRITE(16,24)
      GO TO 2303
      END
      SUBROUTINE CHISQ
      include 'param.inc'
      COMMON/G2/S1,S2,S3,D1,D2,R1,R2,R3,DWD,WDWD
      COMMON/DATAX/Y(IDSZ),YC(IDSZ),KR(IDSZ),BK(IDSZ),NPTS,WT(IDSZ)
      COMMON/JNK/ALOW(30),AHIGH(30),POS(30),BCK(30),
     *NATOM(NUPH),PREF(NUPH,3),NBCKGD,NEXCRG,NSCAT,IWT
      S1=0.0
      S2=0.0
      D1=0.0
      D2=0.0
      S3=0.0
      DWDN=0.0
      DWDD=0.0
      WDWDN=0.0
      WDWDD=0.0
      DO 10 I=1,NPTS
      IF(KR(I).EQ.IRS) GOTO 10
      IF(NBCKGD.NE.0.AND.KR(I).EQ.0) GOTO 10
      DEL=Y(I)-BK(I)-YC(I)
      S1=S1+ABS(DEL)
      S2=S2+DEL*DEL*WT(I)
      D1=D1+Y(I)
      D2=D2+Y(I)*Y(I)*WT(I)
      S3=S3+BK(I)+YC(I)
      IF(I.EQ.1) GO TO 11
      DELK=Y(I-1)-BK(I-1)-YC(I-1)
      DWDEL=( DEL - DELK )**2
      WDWDEL=( DEL * SQRT(WT(I))
     *         - DELK * SQRT(WT(I-1)) )**2
      DWDN=DWDN+DWDEL
      WDWDN=WDWDN+WDWDEL
   11 CONTINUE
      DWDD=DWDD + DEL*DEL
      WDWDD=WDWDD + DEL*DEL*WT(I)
   10 CONTINUE
      R1=0
      R2=100.*S1/D1
      R3=100.*SQRT(S2/D2)
      DWD=DWDN/DWDD
      WDWD=WDWDN/WDWDD
      RETURN
      END
      SUBROUTINE CALCUL(NN)
      include 'param.inc'
      REAL LAMDA
      INTEGER PTR
      COMMON/REFLS/IREFS(IRS),REFS(IRS,8),FF(IRS),IC(NUPH)
      COMMON/G1/THMIN,STEP,THMAX,U,V,W,LAMDA(2),TMV,CTHM,RLIM,SBX,WDT,
     *UU,VV,WW,SYNPOL
      COMMON/CNTRLS/JOBTYP,NPROF,NPHASE,IOT,IPL,IPC,MAT,NXT,LST1,LST2,
     *LST3,IPL1,IPL2,MCYCLE,EPS,MAXS,IFOUR,NFMT
      COMMON/PARAMS/GLB(16),AGLB(16),LGLB(16),PAR(NUPH,20),
     *APAR(NUPH,20),LPAR(NUPH,20),XL(NATS,11),LP(NATS,11),A(NATS,11),
     *PTR(NATS),RELAX(4),RATIO(2)
      CHARACTER*4 ATEXT(NATS),NTYP(NATS)
      COMMON/PARAC/ATEXT,NTYP
      COMMON/JNK/ALOW(30),AHIGH(30),POS(30),BCK(30),
     *NATOM(NUPH),PREF(NUPH,3),NBCKGD,NEXCRG,NSCAT,IWT
      CHARACTER SYMB(NUPH,20)*1
      CHARACTER TITLE*80,PHSNM(NUPH)*80
      COMMON/CHAR/SYMB,TITLE,PHSNM
      COMMON/COEFF/AC(10,16),POSI(30),SCAT(30),DFP(16),DFPP(16)
      CHARACTER*4 NAM(16)
      COMMON/COEFC/NAM
      COMMON/DATAX/Y(IDSZ),YC(IDSZ),KR(IDSZ),BK(IDSZ),NPTS,WT(IDSZ)
      COMMON/RTMTX/IVEC(NUPH,192),MLT(NUPH),ICET(NUPH)
      COMMON/G4/TANN(NIPOI),DERSTO(NIPOI,MSZ)
      COMMON/FOUR/PHASE(IRS)
      DIMENSION HNN(3),XI(14)
      DIMENSION DERIV(MSZ),SNX(NATS),SINTL(NIPOI),SA(NATS),SB(NATS)
      DIMENSION TEMP(NATS),SNEX(NATS),T(3),SUMAX(NATS,9),SUMBX(NATS,9)
      DIMENSION H(3),AL(3,3),SM(3,3)
C     LOGICAL SKIP,PAC
      LOGICAL PAC
      PI=360./3.14159265359
      IPH=IREFS(NN)/256/256/256/4
      IOF=0
      IF(IPH.GT.1) THEN
      DO 101 KK=1,IPH-1
  101 IOF=IOF+NATOM(KK)
      ENDIF
      IRL=MLT(IPH)
      N=NATOM(IPH)
      ICENT=ICET(IPH)
C-----ZEROIZE THE DERIVATIVES OF THIS REFLECTION W.R.T. TO PARAMETERS
      NX=0
      DO 100 I=1,NPHASE
  100 NX=NX+NATOM(I)
      DO 12 I=1,MAXS
   12 DERIV(I)=0.
      DO 21 I=1,NX
      DO 21 J=1,9
      SUMBX(I,J)=0.
   21 SUMAX(I,J)=0.
      CV=0.
      DV=0.
      AV=0.
      BV=0.
      PAC=PAR(IPH,12).NE.0.0.OR.LPAR(IPH,12).NE.0
     *.OR.PAR(IPH,13).NE.0.0.OR.LPAR(IPH,13).NE.0
     *.OR.PAR(IPH,16).NE.0.0.OR.LPAR(IPH,16).NE.0
      DO 303 I=1,3
 303  AL(I,I)=PAR(IPH,I+5)
      AL(1,2)=PAR(IPH,11)
      AL(1,3)=PAR(IPH,10)
      AL(2,1)=0.0
      AL(2,3)=PAR(IPH,9)
      AL(3,1)=0.0
      AL(3,2)=0.0
      AH=PAR(IPH,3)
      BH=PAR(IPH,4)
      CH=PAR(IPH,5)
      B1=PAR(IPH,12)
      B2=PAR(IPH,13)
      UU=PAR(IPH,15)
      VV=PAR(IPH,16)
      WW=PAR(IPH,17)
      NM=MOD(NN,NIPOI)+1
      ICX=MOD(IREFS(NN)/256/256/256,4)
      SLABDA=LAMDA(ICX)*LAMDA(ICX)/4.
      FF(NN)=0.
      PRECOR=1.
      TR=0.
      HNN(3)=MOD(IREFS(NN),256)-128
      HNN(2)=MOD(IREFS(NN)/256,256)-128
      HNN(1)=MOD(IREFS(NN)/256/256,256)-128
      DO 1 I=1,3
    1 H(I)=HNN(I)
C-----CALCULATION OF TEMP.FACTOR,POSITION AND HALFWIDTH
      SS=0.
      DO 2 I=1,3
      DO 2 J=I,3
    2 SS=AL(I,J)*HNN(I)*HNN(J)+SS
      IF(.NOT.PAC)GOTO 13
      TT=0.
      DO 16 I=1,3
      DO 16 J=I,3
   16 TT=TT+AL(I,J)*PREF(IPH,I)*PREF(IPH,J)
      PAK=0.
      DO 33 I=1,3
      DO 33 J=I,3
   33 PAK=AL(I,J)*(PREF(IPH,I)*HNN(J)+PREF(IPH,J)*HNN(I))/2.+PAK
      PAK=PAK*PAK/(TT*SS)
C-----MARCH MODEL FOR PREFERRED ORIENTATION
      PRECOR=(B1*B1*PAK + (1-PAK)/B1)**(-1.5)
   13 SINTL(NM)=SQRT(SS)
      SSNN=0.25*SS
      TAV=EXP(-2.*PAR(IPH,2)*SSNN)*PRECOR
      SINTH=SLABDA*SS
      COSTH=1.-SINTH
      TANTH=SQRT(SINTH/COSTH)
      TANN(NM)=TANTH
      REFS(NN,2)=ATAN(TANTH)*PI
      IF(NPROF.EQ.1.OR.NPROF.EQ.4) GO TO 23
      IF(NPROF.EQ.2) GO TO 24
      REFS(NN,4) = PAR(IPH,15) + SQRT(PAK)*PAR(IPH,16)
      GO TO 20
   23 REFS(NN,4)=WW*REFS(NN,2)*REFS(NN,2)+VV*REFS(NN,2)+UU
      GO TO 20
   24 REFS(NN,4)=WW*REFS(NN,2)*REFS(NN,2)+VV*REFS(NN,2)+UU
      IF(REFS(NN,4).LT.0.6) REFS(NN,4)=0.6
   20 CONTINUE
      CAGLI=AH*TANTH*TANTH+BH*TANTH+CH
C-----FWHM SET TO 0.032 DEGREES WHEN CALCULATED NEGATIVE
      IF (CAGLI.LT.0.0) CAGLI=0.001
      TL=SQRT(CAGLI)
      REFS(NN,1)=TL
      BB=TL*TL
C-----CALCULATION OF COS(H.X),SIN(H.X) AND TEMP. FACTOR FOR EACH ATOM
    3 DO 8 I=1,N
      SNXI=0.
      SAI=0.
      SBI=0.
      DO 9 J=1,11
    9 XI(J)=XL(I+IOF,J)
      DO 10 IR=1,IRL
      DO 301 J=1,3
      IV=IVEC(IPH,IR)/32768/32**(3-J)
      IV=MOD(IV,32)
      SM(J,1)=IV/9-1
      SM(J,2)=MOD(IV/3,3)-1
      SM(J,3)=MOD(IV,3)-1
  301 T(J)=FLOAT(MOD(IVEC(IPH,IR)/32**(3-J),32)-16)/12.0
      X=0.
      DO 6 II=1,3
      YY=0.
      X=T(II)*HNN(II)+X
      DO 7 J=1,3
    7 YY=SM(J,II)*HNN(J)+YY
    6 H(II)=YY
      TR=X
      ARG=TR
      DO 11 J=1,3
   11 ARG=H(J)*XI(J)+ARG
      ARG=6.28318530718*ARG
       ARG2=H(1)*H(1)*XI(6)+H(2)*H(2)*XI(7)+
     1H(3)*H(3)*XI(8)+2.*H(1)*H(2)*XI(9)+2.*
     2 H(1)*H(3)*XI(10)+2.*H(2)*H(3)*XI(11)
      EXPARG=EXP(-ARG2)
      COSA=COS(ARG)*EXPARG
      SINA=SIN(ARG)*EXPARG
      SAI=SAI+COSA
      IF(ICENT.EQ.1)SBI=SINA+SBI
      DO 14 JJ=1,3
      SUMAX(I,JJ)=SUMAX(I,JJ)+H(JJ)*SINA
  14  IF(ICENT.EQ.1)SUMBX(I,JJ)=SUMBX(I,JJ)+H(JJ)*COSA
      SUMAX(I,4)=SUMAX(I,4)+H(1)*H(1)*COSA
      SUMAX(I,5)=SUMAX(I,5)+H(2)*H(2)*COSA
      SUMAX(I,6)=SUMAX(I,6)+H(3)*H(3)*COSA
      SUMAX(I,7)=SUMAX(I,7)+H(1)*H(2)*COSA
      SUMAX(I,8)=SUMAX(I,8)+H(1)*H(3)*COSA
      SUMAX(I,9)=SUMAX(I,9)+H(2)*H(3)*COSA
      IF(ICENT.NE.1)GOTO 10
      SUMBX(I,4)=SUMBX(I,4)+H(1)*H(1)*SINA
      SUMBX(I,5)=SUMBX(I,5)+H(2)*H(2)*SINA
      SUMBX(I,6)=SUMBX(I,6)+H(3)*H(3)*SINA
      SUMBX(I,7)=SUMBX(I,7)+H(1)*H(2)*SINA
      SUMBX(I,8)=SUMBX(I,8)+H(1)*H(3)*SINA
      SUMBX(I,9)=SUMBX(I,9)+H(2)*H(3)*SINA
   10 CONTINUE
       TEMP(I)=EXP(-XL(I+IOF,4)*SSNN)
      SA(I)=SAI
      SB(I)=SBI
      NI=PTR(I+IOF)
      BNI=DFPP(NI)
      FFX=DFP(NI)
      DO 302 II=1,9,2
 302  FFX=FFX+AC(II,NI)*EXP(-AC(II+1,NI)*SSNN)
      SNEXI=FFX*XL(I+IOF,5)*TEMP(I)
      SNXI=BNI*XL(I+IOF,5)*TEMP(I)
C-----CALCULATE A AND B OF F
      AV=SNEXI*SAI+AV
      BV=SNEXI*SBI+BV
      SNEX(I)=2.*SNEXI*TAV*RATIO(ICX)
      SNX(I)=2.*SNXI*TAV*RATIO(ICX)
      CV=CV+SNXI*SAI
      DV=DV+SNXI*SBI
    8 CONTINUE
      FNN=RATIO(ICX)*(CV*CV+AV*AV+DV*DV+BV*BV)*TAV
  120 FF(NN)=FNN
      PHASE(NN)=1.57079633
      IF(AV.NE.0.0) PHASE(NN)=ATAN(BV/AV)
      IF(AV.LT.0.0) PHASE(NN)=PHASE(NN)+3.14159265
      IF(PHASE(NN).LT.0.0) PHASE(NN)=PHASE(NN)+6.28318530
      IF(MAXS.EQ.0)RETURN
C-----CALCULATE DERIVATIVES
      DO 17 I=1,N
      SNEXI=SNEX(I)
      SNXI=SNX(I)
      SAI=SA(I)
      SBI=SB(I)
      DO 22 J=1,11
      K=LP(I+IOF,J)
      IF(K.EQ.0) GOTO22
      IF(J.GT.5) GOTO 221
      IF(J.GT.3)GOTO 29
      SUMA=SUMAX(I,J)
      SUMB=SUMBX(I,J)
      DER=-((AV*SUMA-BV*SUMB)*SNEXI+(CV*SUMA-DV*SUMB)*SNXI)
     1      *6.28318530718
      GOTO 26
   29 IF(J.GT.4)GOTO 31
      DER=-((SAI*AV+SBI*BV)*SNEXI+(SAI*CV+SBI*DV)*SNXI)*SSNN
      GOTO 26
   31 DER=((SAI*AV+SBI*BV)*SNEXI+(SAI*CV+SBI*DV)*SNXI)/XL(I+IOF,5)
      GOTO 26
  221 SUMA=SUMAX(I,J-2)
      SUMB=SUMBX(I,J-2)
      DER=-((AV*SUMA+BV*SUMB)*SNEXI+(CV*SUMA+DV*SUMB)*SNXI)
      IF(J.GE.9)DER=2.*DER
   26 DERIV(K)=A(I+IOF,J)*DER+DERIV(K)
   22 CONTINUE
   17 CONTINUE
C-----CALCULATE DERIVATIVES
C-----MARCH MODEL FOR PREFERRED ORIENTATION
      K=LPAR(IPH,12)
      IF(K.NE.0)DERIV(K)=DERIV(K)-APAR(IPH,12)*1.5*(B1*B1*PAK + (1-PAK
     & )/B1)**(-2.5)*(2*B1*PAK - (1-PAK)/B1**2) *FNN/PRECOR
      K=LPAR(IPH,2)
      IF(K.NE.0)DERIV(K)=DERIV(K)-2.*APAR(IPH,2)*SSNN*FNN
      K=LPAR(IPH,1)
      IF(K.NE.0)DERIV(K)=DERIV(K)+APAR(IPH,1)*FNN/PAR(IPH,1)
C-----DERIVATIVES FOR VARIABLE SHAPE IN PSEUDO-VOIGT AND PEARSON VII
      X=1.0
      DO 15 J=15,17
      K=LPAR(IPH,J)
      IF(K.EQ.0) GO TO 19
      IF(NPROF.NE.3.AND.NPROF.NE.5) GO TO 18
      IF(J.EQ.15) DERIV(K)=DERIV(K)+APAR(IPH,J)*FNN
      IF(J.EQ.16) DERIV(K)=DERIV(K)+SQRT(PAK)*APAR(IPH,J)*FNN
      GO TO 15
   18 DERIV(K)=DERIV(K)+X*APAR(IPH,J)*FNN
   19 X=X*REFS(NN,2)
   15 CONTINUE
C     END CODE
      SINTH=0.5*FNN*PI*SLABDA/SQRT(SINTH*COSTH)
      SS=0.5*FNN/TL
      X=TANTH*TANTH
      DO 78 J=3,5
      K=LPAR(IPH,J)
      IF(K.EQ.0)GOTO 78
      DERIV(K)=X*APAR(IPH,J)*SS+DERIV(K)
   78 X=X/TANTH
      K=LGLB(1)
      IF(K.NE.0) DERIV(K)=AGLB(1)*FNN+DERIV(K)
      DO 79 J=1,6
      K=LPAR(IPH,J+5)
      IF(K.EQ.0)GOTO 79
      IF(J.LT.4)X=HNN(J)*HNN(J)
      IF(J.EQ.4)X=HNN(2)*HNN(3)
      IF(J.EQ.5)X=HNN(1)*HNN(3)
      IF(J.EQ.6)X=HNN(1)*HNN(2)
      DERIV(K)=X*SINTH*APAR(IPH,J+5)+DERIV(K)
   79 CONTINUE
      K=LPAR(IPH,14)
      IF(K.NE.0) DERIV(K)=APAR(IPH,14)*FNN+DERIV(K)
C-----STORE DERIVATIVES FOR LIMO REFLECTIONS AT A TIME
      DO 80 I=1,MAXS
   80 DERSTO(NM,I)=DERIV(I)
   48 RETURN
      END
      SUBROUTINE REWRIT(ISCALE,IDIF,IFY,IFYC,IFM,IFD,EXPAND)
      include 'param.inc'
      REAL LAMDA
      INTEGER PTR
      COMMON/DC/SAVE(NUPH,8),NSAVE
      COMMON/G1/THMIN,STEP,THMAX,U,V,W,LAMDA(2),TMV,CTHM,RLIM,SB,WDT,
     *UU,VV,WW,SYNPOL
      COMMON/CNTRLS/JOBTYP,NPROF,NPHASE,IOT,IPL,IPC,MAT,NXT,LST1,LST2,
     *LST3,IPL1,IPL2,MCYCLE,EPS,MAXS,IFOUR,NFMT
      COMMON/PARAMS/GLB(16),AGLB(16),LGLB(16),PAR(NUPH,20),
     *APAR(NUPH,20),LPAR(NUPH,20),XL(NATS,11),LP(NATS,11),A(NATS,11),
     *PTR(NATS),RELAX(4),RATIO(2)
      CHARACTER*4 ATEXT(NATS),NTYP(NATS)
      COMMON/PARAC/ATEXT,NTYP
      COMMON/JNK/ALOW(30),AHIGH(30),POS(30),BCK(30),
     *NATOM(NUPH),PREF(NUPH,3),NBCKGD,NEXCRG,NSCAT,IWT
      CHARACTER SYMB(NUPH,20)*1
      CHARACTER TITLE*80,PHSNM(NUPH)*80
      COMMON/CHAR/SYMB,TITLE,PHSNM
      COMMON/COEFF/AC(10,16),POSI(30),SCAT(30),DFP(16),DFPP(16)
      CHARACTER*4 NAM(16)
      COMMON/COEFC/NAM
c     OPEN(UNIT=8,STATUS='NEW',CARRIAGE CONTROL='LIST')
      WRITE(8,1)TITLE
   1  FORMAT(A80)
      JOBTYP=JOBTYP-1
      WRITE(8,2)JOBTYP,NPROF,NPHASE,NBCKGD,NEXCRG,NSCAT,IWT
   2  FORMAT(9I4)
      WRITE(8,14)IOT,IPL,IPC,MAT,NXT,LST1,LST2,LST3,IPL1,IPL2
     *,IFOUR,NFMT
  14  FORMAT(12I1)
      WRITE(8,23)LAMDA,RATIO(2),SB,WDT,CTHM,TMV,RLIM,SYNPOL
  23  FORMAT(2F8.5,7F8.4)
      WRITE(8,30)MCYCLE,EPS,RELAX,THMIN,STEP,THMAX
  30  FORMAT(1I4,5F4.2,3F8.3)
      IF(NBCKGD.LT.2)GOTO 120
      WRITE(8,34)(POS(I),BCK(I),I=1,NBCKGD)
  34  FORMAT(2F8.2)
 120  IF(NEXCRG.LE.0)GOTO 122
      WRITE(8,34)(ALOW(I),AHIGH(I),I=1,NEXCRG)
 122  IF(NSCAT.LE.0)GOTO 124
      DO 125 I=1,NSCAT
      WRITE(8,38)NAM(I),DFP(I),DFPP(I)
 38   FORMAT(A4,2F8.4)
 126  IF(JOBTYP.EQ.0)WRITE(8,39)(AC(J,I),J=1,9)
 39   FORMAT(F8.5,8F9.5)
 125  CONTINUE
  124 CONTINUE
      WRITE(8,43)MAXS
   43 FORMAT(I8)
      N=0
      DO 121 I=1,NPHASE
  121 N=N+NATOM(I)
      DO 84 I=1,N
      DO 84 J=1,11
  84  A(I,J)=SIGN(1.,A(I,J))*(FLOAT(10*LP(I,J))+ABS(A(I,J)))
      DO 87 I=1,NPHASE
      DO 88 J=1,6
   88 PAR(I,J+5)=SAVE(I,J)
      DO 87 J=1,17
  87  APAR(I,J)=SIGN(1.,APAR(I,J))*(FLOAT(10*LPAR(I,J))+ABS(APAR(I,J)))
      DO 89 J=1,7
  89  AGLB(J)=SIGN(1.,AGLB(J))*(FLOAT(10*LGLB(J))+ABS(AGLB(J)))
      WRITE(8,45)GLB(1),AGLB(1)
   45 FORMAT(2F8.4)
      IF(NBCKGD.EQ.0)WRITE(8,46)(GLB(J),J=2,7),(AGLB(J),J=2,7)
   46 FORMAT(6G12.6,/,6F12.6)
 477  DO 81 K=1,NPHASE
      IOF=0
      IF(K.GT.1) THEN
      DO 100 KK=1,K-1
  100 IOF=IOF+NATOM(KK)
      ENDIF
      WRITE(8,465)PHSNM(K)
  465 FORMAT(A80)
      WRITE(8,48)NATOM(K),(PREF(K,I),I=1,3)
   48 FORMAT(I4,8X,3F4.1)
      N=NATOM(K)
      WRITE(8,460)(SYMB(K,I),I=1,20)
  460 FORMAT(20A1)
      WRITE(8,65)(ATEXT(I+IOF),NTYP(I+IOF),(XL(I+IOF,J),J=1,11),I=1,N)
   65 FORMAT(2A4,8X,5F8.5/6F8.5)
      WRITE(8,50)(PAR(K,I),I=1,17)
   50 FORMAT(G12.6,4X,F8.5,/,3F8.5,/,6F8.4,/,3F8.5,/,F8.5,F8.5,F8.6)
      WRITE(8,861)((A(I+IOF,J),J=1,11),I=1,N)
  861 FORMAT(5F8.3,/,6F8.3)
      WRITE(8,86)(APAR(K,I),I=1,17)
   86 FORMAT(F12.6,4X,F8.3,/,3F8.3,/,6F8.3,/,3F8.3,/,3F8.3)
   81 CONTINUE
      IF(IPL.NE.0)WRITE(8,191)ISCALE,IDIF
  191 FORMAT(2I8)
      IF(IPL2.NE.0)WRITE(8,171)IFY,IFYC,IFM,IFD,EXPAND
  171 FORMAT(4I1,4X,F8.4)
  151 RETURN
      END
      SUBROUTINE OUTPTR(ICYCLE)
      include 'param.inc'
      REAL LAMDA
      COMMON/DC/SAVE(NUPH,8),NSAVE
      COMMON/PARAMS/GLB(16),AGLB(16),LGLB(16),PAR(NUPH,20),
     *APAR(NUPH,20),LPAR(NUPH,20),XL(NATS,11),LP(NATS,11),A(NATS,11),
     *PTR(NATS),RELAX(4),RATIO(2)
      CHARACTER*4 ATEXT(NATS),NTYP(NATS)
      COMMON/PARAC/ATEXT,NTYP
      COMMON/JNK/ALOW(30),AHIGH(30),POS(30),BCK(30),
     *NATOM(NUPH),PREF(NUPH,3),NBCKGD,NEXCRG,NSCAT,IWT
      CHARACTER SYMB(NUPH,20)*1
      CHARACTER TITLE*80,PHSNM(NUPH)*80
      COMMON/CHAR/SYMB,TITLE,PHSNM
      COMMON/G3/COND,IORD1,IORD2,TH,NUM
      COMMON/CNTRLS/JOBTYP,NPROF,NPHASE,IOT,IPL,IPC,MAT,NXT,LST1,LST2,
     *LST3,IPL1,IPL2,MCYCLE,EPS,MAXS,IFOUR,NFMT
      COMMON/G2/S1,S2,S3,D1,D2,R1,R2,R3,DWD,WDWD
      COMMON/G1/THMIN,STEP,THMAX,FILL(3),LAMDA(2),TMV,CTHM,RLIM,SB,WDT,
     *UU,VV,WW,SYNPOL
      COMMON/F1/SMM(MSZ,MSZ),V(MSZ)
      DIMENSION SY(20),SZ(20),DCSM(6,6),DCV(6)
C***
      WRITE(16,1) ICYCLE
      WRITE(7,111) ICYCLE
 111  FORMAT(//14H CYCLE NUMBER=,I4)
    1 FORMAT(1H ,120(1H+)//14H CYCLE NUMBER=,I4)
      DO 10 IP=1,NPHASE
      IOF=0
      IF(IP.NE.1) THEN
      DO 104 KK=1,IP-1
  104 IOF=IOF+NATOM(KK)
      ENDIF
      WRITE(16,3)IP,PHSNM(IP)
    3 FORMAT(7H PHASE ,I2,2H: ,A80,/' NEW PARAMETERS, SHIFTS, AND ',
     *'STANDARD DEVIATIONS='//5H ATOM,6X,1HX,7X,2HDX,6X,2HSX,7X,1HY,7X,
     *2HDY,6X,2HSY,7X,1HZ,7X,2HDZ,6X,2HSZ,6X,1HB,6X,2HDB,5X,2HSB,
     *6X,1HN,6X,2HDN,5X,2HSN)
      N=NATOM(IP)
      DO  4 I=1,N
      DO  5 J=1,5
      KM=LP(I+IOF,J)
      IF(KM.NE.0) GOTO  6
      SY(J)=0.
      SZ(J)=0.
      GOTO  5
    6 SZ(J)=SQRT(ABS(SMM(KM,KM)))*ABS(A(I+IOF,J))
      SY(J)=V(KM)*A(I+IOF,J)*RELAX(1)
      XL(I+IOF,J)=XL(I+IOF,J)+SY(J)
    5 CONTINUE
    4 WRITE(16,100)ATEXT(I+IOF),(XL(I+IOF,J),SY(J),SZ(J),J=1,5)
  100 FORMAT(1X,A4,2X,3(3F8.5,1X),2(3F7.4,1X))
      WRITE(16,101)
  101 FORMAT(5H ATOM,7X,3HB11,6X,4HDB11,6X,4HSB11,7X,3HB22,6X,4HDB22,
     !6X,4HSB22,7X,3HB33,6X,4HDB33,6X,4HSB33,/,12X,3HB12,6X,4HDB12,6X,
     !4HSB12,7X,3HB13,6X,4HDB13,6X,4HSB13,7X,3HB23,6X,4HDB23,6X,4HSB23)
      DO  7 I=1,N
      DO  8 J=6,11
      KM=LP(I+IOF,J)
      IF(KM.NE.0)GOTO  9
      SY(J)=0.
      SZ(J)=0.
      GOTO  8
    9 SZ(J)=SQRT(ABS(SMM(KM,KM)))*ABS(A(I+IOF,J))
      SY(J)=V(KM)*A(I+IOF,J)*RELAX(2)
      XL(I+IOF,J)=XL(I+IOF,J)+SY(J)
    8 CONTINUE
    7 WRITE(16,102)ATEXT(I+IOF),(XL(I+IOF,J),SY(J),SZ(J),J=6,11)
  102 FORMAT(1X,A4,9F10.6,/,5X,9F10.6)
      DO 17 J=1,17
      KM=LPAR(IP,J)
      IF(KM.NE.0)GOTO 18
      SY(J)=0.
      SZ(J)=0.
      GOTO 17
   18 SZ(J)=SQRT(ABS(SMM(KM,KM)))*ABS(APAR(IP,J))
      SY(J)=V(KM)*APAR(IP,J)*RELAX(3)
      PAR(IP,J)=PAR(IP,J)+SY(J)
   17 CONTINUE
      CALL DIRECT(DCSM,DCV,IP,N)
      DO 30 I=1,6
      SY(I+5)=DCV(I)-SAVE(IP,I)
   30 SAVE(IP,I)=DCV(I)
      IF(NPROF.EQ.3.OR.NPROF.EQ.5) GO TO 200
      WRITE(16,103)(PAR(IP,I),SY(I),SZ(I),I=1,5),
     *(DCV(I),SY(I+5),DCSM(I,I),I=1,6),(PAR(IP,I),SY(I),SZ(I),
     *I=12,17)
  103 FORMAT(' OVERALL SCALE FACTOR=',1X,G12.6,1X,G9.3,1X,G9.3,/,
     *' OVERALL TEMP. FACTOR=',3F10.4,/,
     *' HALFWIDTH PARAMETERS=',/,3(1X,3F10.6,/),
     *' CELL PARAMETERS=',/,6(3F11.6,/),
     *' PREFERRED ORIENTATION PARAMETERS=',/,2(3F9.5,/),
     *' ASYMMETRY PARAMETER=',/,3F8.5,/,
     *' PEAK SHAPE PARAMETERS = ',/,3F10.4,/,3F10.5,/,3F10.6,/)
      GO TO 201
  200 WRITE(16,202)(PAR(IP,I),SY(I),SZ(I),I=1,2),
     *(DCV(I),SY(I+5),DCSM(I,I),I=1,6),
     *(PAR(IP,I),SY(I),SZ(I),I=12,14),
     *(PAR(IP,I),SY(I),SZ(I),I=3,5),PAR(IP,15),SY(15),SZ(15),
     *PAR(IP,16),SY(16),SZ(16)
  202 FORMAT(' OVERALL SCALE FACTOR = ',1X,G12.6,1X,G9.3,1X,G9.3,/,
     *' OVERALL TEMP. FACTOR = ',3F10.4,/,
     *' CELL PARAMETERS = ',/,6(3F11.6,/),
     *' PREFERRED ORIENTATION PARAMETERS = ',/,2(3F9.5,/),
     *' ASYMMETRY PARAMETER = ',/,1X,3F8.5,/,
     *' VOIGT PEAK WIDTH/SHAPE :',
     */,1X,'GAUSSIAN COMPONENT HALF-WIDTH PARAMETERS = ',/,
     *3(1X,3F10.6,/),1X,'LORENTZIAN COMPONENT PARAMETERS = ',/,
     *2(1X,3F10.6,/))
      IF(PAR(IP,15).NE.0.) THEN
      PARSIZ=57.29577951*LAMDA(1)/PAR(IP,15)
      WRITE(16,70) PARSIZ
      ENDIF
      IF(PAR(IP,16).NE.0.) THEN
      PARSIZ=57.29577951*LAMDA(1)/(PAR(IP,15)+PAR(IP,16))
      WRITE(16,71) PARSIZ
      ENDIF
   70 FORMAT(' EQUIVALENT TO A PARTICLE SIZE OF',F10.1,' ANGSTROMS')
   71 FORMAT('                               BY',F10.1,' ANGSTROMS')
  201 CONTINUE
   10 CONTINUE
      DO 19 J=1,7
      KM=LGLB(J)
      IF(KM.NE.0)GOTO 20
      SY(J)=0.
      SZ(J)=0.
      GOTO 19
   20 SZ(J)=SQRT(ABS(SMM(KM,KM)))*ABS(AGLB(J))
      SY(J)=V(KM)*AGLB(J)*RELAX(4)
      GLB(J)=GLB(J)+SY(J)
   19 CONTINUE
      WRITE(16,109)GLB(1),SY(1),SZ(1)
  109 FORMAT(' GLOBAL PARAMETERS',/,' ZEROPOINT=',3F8.4)
      IF(NBCKGD.NE.0)GOTO 80
      WRITE(16,84)(GLB(J),SY(J),SZ(J),J=2,7)
   84 FORMAT(' BACKGROUND PARAMETERS',/,(1X,3G12.6))
   80 CONTINUE
      WRITE(16,57)R2,R3
      X=100.*SQRT((FLOAT(NUM)-FLOAT(MAXS))*1./D2)
      WRITE(16,58)X
      WRITE(16,64)DWD,WDWD
   64 FORMAT(29H DURBIN-WATSON D STATISTIC : ,/,
     + 10X,'UNWEIGHTED = ',F8.3,/,10X,'WEIGHTED =   ',F8.3)
   60 I=NUM-MAXS
      WRITE(16,62)I
      S2=S2/FLOAT(I)
      WRITE(16,63)S1,D1,S3,D2,S2,COND
      WRITE(7,59) R2,R3,S2
   59 FORMAT(' RP = ',F6.2,'   RWP = ',F6.2,'   GOF = ',F15.2)
   57 FORMAT( 12H R-FACTORS =,2F8.2)
   58 FORMAT(11H EXPECTED =,F8.2)
   62 FORMAT( 7H N-P  =,I8)
   63 FORMAT(/,7X,7HSUMYDIF, 7X,7HSUMYOBS, 6X,8HSUMYCALC, 6X,10HSUMWYOBS
     1SQ,6X,8HG OF FIT,6X,9HCONDITION/1X,6E14.4/)
      RETURN
      END
