                    SUBROUTINE MEDWT(MEMSZ,DELMAX)                              
C                                                                               
C       THIS SUBROUTINE DETERMINES THE MEDIAN VALUE OF THE FUNCTION             
C             ABS (FOBS,FCALC,SIGMA)=ABS((FOBS-FCALC)/SIGMA)                    
C       FROM THE PREVIOUS REFINEMENT CYCLE. THIS VALUE IS MULTIPLIED            
C              9 AND THE RESULT IS PLACED IN DELMAX.                            
C                                                                               
C                      WRITTEN BY L. W. FINGER                                  
C                                                                               
C                                                                               
      DIMENSION RMAT(3,3,24),TRANS(3,24),IFBUF(504)                             
C                                                                               
      COMMON /LSMAT/ B(200),IEXT,IPARA(1000),MRANK,KCYL,LABEL(3),JSCRF          
     & (10), NRSF, SIGMA,MDEGN,MAXIMA,IALFA(60),RDAMP(200),IDAMP(200),          
     & A(2200)                                                                  
      COMMON /TAPE/ IN,INN,IOUT,ISC1,ISC2,ISC3,ISC4,ISC5,IPUN                   
      COMMON /FBUF/ FBUF(504)                                                   
C                                                                               
      EQUIVALENCE (FBUF(1),IFBUF(1),VOL,RMAT(1,1,1),TRANS(1,1))                 
C                                                                               
    2 FORMAT(24H0MEDIAN DELTA F/SIGMA IS,F9.4,14H NEW DELMAX IS,F9.4)           
   93 FORMAT(13H0ERRORMESSAGE/1X,12(1H-)//                                      
     &62H SUBROUTINE MEDWT HAS TRIED TO USE THE ARRAY A (COMMON LSMAT),/        
     &62H0BUT THE DIMENSION IS INSUFFICIENT.        FOR CALCULATION OF /        
     &62H0CORRECT SIZE SEE GENERAL DESCRIPTION OF PROMETHEUS.          /        
     &16H0JOB TERMINATED.)                                                      
      REWIND ISC4                                                               
      READ(ISC4)NSYM,VOL,(((RMAT(I,J,K),I=1,3),J=1,3),K=1,NSYM)                 
     & ,((TRANS(I,J),I=1,3),J=1,NSYM)                                           
      NREF=0                                                                    
   10 IF(MEMSZ-MDEGN-NREF-50)11,15,15                                           
   11 MEMSZ=MEMSZ+256                                                           
      IF(MAXIMA-MEMSZ)91,92,92                                                  
   91 WRITE(IOUT,93)                                                            
      CALL EXIT                                                                 
C                                                                               
   92 CONTINUE                                                                  
   15 READ(ISC4)FBUF                                                            
      M=1                                                                       
   20 IF(IFBUF(M).EQ.99)GO TO 100                                               
      IF(IFBUF(M+7).NE.0)GO TO 25                                               
      NREF=NREF+1                                                               
      A(MDEGN+NREF)=ABS((FBUF(M+3)-FBUF(M+4))/FBUF(M+8))                        
   25 M=M+10                                                                    
      IF(M-500)20,20,10                                                         
C---- ALL REFLECTIONS READ, NOW SORT THEM                                       
  100 IX=MDEGN+1                                                                
      CALL FNDMED(A(IX),NREF,XMED)                                              
  150 DELMAX=9.0*XMED                                                           
      WRITE(IOUT,2)XMED,DELMAX                                                  
      RETURN                                                                    
C                                                                               
      END                                                                       
C                                                                               
                     SUBROUTINE FNDMED(X,N,XMED)                                
C                                                                               
C     THIS SUBROUTINE FINDS THE MEDIAN OF THE N ELEMENTS OF THE ARRAY           
C       IF N IS EVEN THE DEDIAN IS TAKEN TO BE THE MEAN OF THE TWO              
C                        MIDDLE ELEMENTS.                                       
C                                                                               
C                        BY L. W. FINGER                                        
      DIMENSION X(N)                                                            
C                                                                               
      COMMON /TAPE/ IN,INN,IOUT,ISC1,ISC2,ISC3,ISC4,ISC5,IPUN                   
 9000 FORMAT(13H0ERRORMESSAGE/1X,12(1H-)//                                      
     & 47H NULL LIST HAS NO MEDIAN. (DELMAX COMPUTATION.)/                      
     & 16H0JOB TERMINATED.)                                                     
C                                                                               
      IF(N.GT.0)GO TO 10                                                        
      WRITE(IOUT,9000)                                                          
      CALL EXIT                                                                 
C                                                                               
   10 IF(N.GE.3)GO TO 20                                                        
      IF(N.EQ.1)XMED=X(1)                                                       
      IF(N.EQ.2)XMED=(X(1)+X(2))*0.5                                            
      RETURN                                                                    
C                                                                               
   20 MEDHI=N/2+1                                                               
      MEDLO=MEDHI                                                               
      IF(MOD(N,2).EQ.0)MEDLO=MEDLO-1                                            
C---- DETERMINE IF X IS ALREADY SORTED                                          
      M=N-1                                                                     
      DO 30 I=1,M                                                               
      IF(X(I).GT.X(I+1))GO TO 40                                                
   30 CONTINUE                                                                  
   35 XMED=(X(MEDLO)+X(MEDHI))*0.5                                              
      RETURN                                                                    
C                                                                               
C---- LIST NOT ALREADY SORTED.  IT MUST BE PARTIALLY SORTED TO                  
C     FIND THE MEDIAN.                                                          
   40 I=1                                                                       
      J=N                                                                       
   45 IF(J-I-9)100,100,46                                                       
   46 K=I                                                                       
      MID=(I+J)/2                                                               
      AMED=X(MID)                                                               
      IF(X(I).LE.AMED)GO TO 50                                                  
   47 X(MID)=X(I)                                                               
      X(I)=AMED                                                                 
      AMED=X(MID)                                                               
   50 L=J                                                                       
      IF(X(J).GE.AMED)GO TO 60                                                  
      X(MID)=X(J)                                                               
      X(J)=AMED                                                                 
      AMED=X(MID)                                                               
      IF(X(I).GT.AMED)GO TO 47                                                  
      GO TO 60                                                                  
   55 X(L)=X(K)                                                                 
      X(K)=TT                                                                   
C---- AMED IS A GUESS AT THE MEDIAN OF THE CURRENT SET. DIVIDE IT INTO          
C     TWO SUBSETS THAT ARE GREATER THAN OR LESS THAN AMED.                      
   60 L=L-1                                                                     
      IF(L.LE.K)GO TO 70                                                        
      IF(X(L).GT.AMED)GO TO 60                                                  
      TT=X(L)                                                                   
   65 K=K+1                                                                     
      IF(K.GE.L)GO TO 70                                                        
      IF(X(K).LT.AMED)GO TO 65                                                  
      GO TO 55                                                                  
   70 L=K+1                                                                     
C---- K NOW POINTS AT TOP OF LOW SUBSET, L AT BOTTOM OF HIGH SUBSET.            
   72 IF(K.LE.MEDHI)GO TO 75                                                    
      J=K                                                                       
      GO TO 45                                                                  
   75 IF(L.GE.MEDLO)GO TO 80                                                    
C---- MEDIAN IS IN INTERIOR OF HIGH SUBSET                                      
      I=L                                                                       
      GO TO 45                                                                  
   80 IF(K.NE.MEDLO)GO TO 85                                                    
C---- K=MEDLO,  X(MEDLO)=MAXIMUM OF LOW SUBSET.                                 
      MM=MEDLO                                                                  
   82 AMED=X(MM)                                                                
      NN=MM-1                                                                   
      DO 84 II=I,NN                                                             
      IF(X(II).LT.AMED)GO TO 84                                                 
      X(MM)=X(II)                                                               
      X(II)=AMED                                                                
      AMED=X(MM)                                                                
   84 CONTINUE                                                                  
      K=K-1                                                                     
      GO TO 80                                                                  
   85 IF(K.NE.MEDHI)GO TO 90                                                    
C---- K=MEDHI,  X(MEDHI)=MINIMUM OF SUBSET.                                     
      MM=MEDHI                                                                  
      GO TO 82                                                                  
   90 IF(L.NE.MEDHI)GO TO 95                                                    
C---- L=MEDHI,  X(MEDHI)=MINIMUM OF HIGH SUBSET.                                
      MM=MEDHI                                                                  
   92 AMED=X(MM)                                                                
      NN=MM+1                                                                   
      DO 94 II=NN,J                                                             
      IF(X(II).GT.AMED)GO TO 94                                                 
      X(MM)=X(II)                                                               
      X(II)=AMED                                                                
      AMED=X(MM)                                                                
   94 CONTINUE                                                                  
      L=L+1                                                                     
      GO TO 90                                                                  
   95 IF(L.NE.MEDLO)GO TO 35                                                    
C---- L=MEDLO,  X(MEDLO)=MINIMUM OF HIGH SUBSET                                 
      MM=MEDLO                                                                  
      GO TO 92                                                                  
C---- A SUBSET OF 9 OR FEWER ELEMENTS CONTAINS THE MEDIAN. ORDER THE            
C     SUBSET BY BRUTE FORCE.                                                    
  100 K=I-1                                                                     
  102 K=K+1                                                                     
      IF(K.EQ.J)GO TO 35                                                        
      AMED=X(K+1)                                                               
      IF(X(K).LE.AMED)GO TO 102                                                 
      L=K                                                                       
  105 X(L+1)=X(L)                                                               
      IF(L.EQ.I.OR.AMED.GE.X(L-1))GO TO 110                                     
      L=L-1                                                                     
      GO TO 105                                                                 
  110 X(L)=AMED                                                                 
      GO TO 102                                                                 
      END                                                                       
C                                                                               
                      SUBROUTINE JACOBI(QTQ)                                    
C                                                                               
C     CALCULATION OF EIGENVALUES OF THE ORTHOGONAL EXTINCTION TENSOR.           
C          PROGRAMME ADAPTED FROM XRAY SYSTEM  BY U. ZUCKER                     
C                                                                               
      DIMENSION QTQ(100)                                                        
C                                                                               
      COMMON/TAPE/IN,INN,IOUT,ISC1,ISC2,ISC3,ISC4,ISC5,IPUN                     
C                                                                               
   13 FORMAT(43H0CAUTION - TROUBLES IN DIAGONALIZING TENSOR)                    
C                                                                               
  140 QTQ(19)=0.0                                                               
      I=0                                                                       
      DO 1519 J=1,3                                                             
      DO 1519 K=1,3                                                             
      I=I+1                                                                     
      IF(J-K)1513,1515,1513                                                     
 1513 QTQ(I+9)=0.0                                                              
      QTQ(19)=QTQ(19)+QTQ(I)**2                                                 
      GOTO 1519                                                                 
 1515 QTQ(I+9)=1.0                                                              
      QTQ(19)=QTQ(19)+QTQ(I)                                                    
 1519 CONTINUE                                                                  
      QTQ(19)=QTQ(19)/3.                                                        
      QTQ(20)=0.0001*QTQ(19)                                                    
 1531 QTQ(19)=0.333333*QTQ(19)                                                  
      I10=0                                                                     
 1533 IT=0                                                                      
      I10=I10+1                                                                 
      IF(I10-10)11,11,12                                                        
   12 WRITE(IOUT,13)                                                            
      GOTO 1595                                                                 
   11 CONTINUE                                                                  
      JUNK1=1                                                                   
 1535 JUNK1=JUNK1+1                                                             
      JUNK2=0                                                                   
 1537 JUNK2=JUNK2+1                                                             
      I=3*(JUNK1-1)+JUNK2                                                       
      IF(ABS(QTQ(I))-QTQ(19))1591,1551,1551                                     
 1551 J=3*(JUNK1-1)+JUNK1                                                       
      K=3*(JUNK2-1)+JUNK2                                                       
      L=3*(JUNK2-1)+JUNK1                                                       
      IT=IT+1                                                                   
      QTQ(21)=0.5*(QTQ(K)-QTQ(J))                                               
      TEMP=SQRT(QTQ(L)**2+QTQ(21)**2)                                           
      IF(QTQ(21))1555,1553,1553                                                 
 1553 QTQ(22)=-QTQ(L)/TEMP                                                      
      GOTO 1557                                                                 
 1555 QTQ(22)=QTQ(L)/TEMP                                                       
 1557 QTQ(22)=SIGN(AMIN1(ABS(QTQ(22)),1.0),QTQ(22))                             
      SINR=QTQ(22)/(SQRT(2.0*(1.0+SQRT(1.0-QTQ(22)**2))))                       
      SINSQ=SINR**2                                                             
      COSSQ=1.0-SINSQ                                                           
      COSR=SQRT(COSSQ)                                                          
      SINCOS=SINR*COSR                                                          
      QTQ(23)=QTQ(J)                                                            
      QTQ(24)=QTQ(K)                                                            
      QTQ(25)=QTQ(L)                                                            
      DO 1579 M=1,3                                                             
      I1=3*(JUNK2-1)+M                                                          
      I2=3*(JUNK1-1)+M                                                          
      I3=3*(M-1)+JUNK2                                                          
      I4=3*(M-1)+JUNK1                                                          
      QTQ(26)=QTQ(I1+9)                                                         
      QTQ(I1+9)=QTQ(I1+9)*COSR-QTQ(I2+9)*SINR                                   
      QTQ(I2+9)=QTQ(26)*SINR+QTQ(I2+9)*COSR                                     
      IF(M-JUNK2) 1571,1579,1571                                                
 1571 IF(M-JUNK1) 1573,1579,1573                                                
 1573 QTQ(26)=QTQ(I3)                                                           
      QTQ(I3)=QTQ(I3)*COSR-QTQ(I4)*SINR                                         
      QTQ(I4)=QTQ(26)*SINR+QTQ(I4)*COSR                                         
      QTQ(I1)=QTQ(I3)                                                           
      QTQ(I2)=QTQ(I4)                                                           
 1579 CONTINUE                                                                  
      QTQ(K)=QTQ(24)*COSSQ+QTQ(23)*SINSQ-2.0*QTQ(25)*SINCOS                     
      QTQ(J)=QTQ(24)*SINSQ+QTQ(23)*COSSQ+2.0*QTQ(25)*SINCOS                     
      QTQ(L)=(QTQ(24)-QTQ(23))*SINCOS+QTQ(25)*(COSSQ-SINSQ)                     
      QTQ(I)=QTQ(L)                                                             
 1591 IF(JUNK2-JUNK1+1)1537,1593,1593                                           
 1593 IF(JUNK1-3)1535,1594,1594                                                 
 1594 IF(IT)1595,1595,1533                                                      
 1595 IF(QTQ(19)-QTQ(20))1597,1597,1531                                         
 1597 RETURN                                                                    
C                                                                               
      END                                                                       
                                                                                
                                                                                
                                                                                
                                                                                
                                                                                
                                                                                
                                                                                
