//BERMUDA  JOB                                                                  
/*JOBPARM S=3,R=2048,L=1                                                        
//STEP1   EXEC FORTVCLG,RC=1024K,RL=1024K,RG=3072K,                             
//        PARM.FORT='LANGLVL(77),NOMAP,NOXREF,OPTIMIZE(1)',                     
//        PARM.LKED='LET'                                                       
//FORT.SYSIN DD *                                                               
      IMPLICIT REAL*8(A-H,O-Z)                                                  
      DIMENSION ARRAY(1000), QUANT(1000) ,PROB(1000)                            
      IPRINT = 6                                                                
      IDISK  =10                                                                
      IREAD  = 5                                                                
C                                                                               
C     READ TEST DATA                                                            
C---------------------------------------------------                            
      READ(IREAD,*) NDATA                                                       
      DO 1 I=1,NDATA                                                            
      READ(IREAD,*) ARRAY(I)                                                    
 1    CONTINUE                                                                  
C                                                                               
C-----------------------------------------------------------                    
C     PERFORM Q-Q PLOTS FOR DEPARTURE FROM NORMALITY                            
C-----------------------------------------------------------                    
      CALL QQTEST(ARRAY,NDATA,IPRINT,IDISK,QUANT,PROB,IERR,RQ)                  
C                                                                               
C-----------------------------------------------------------                    
C     WRITE DATA ON DISK                                                        
C-----------------------------------------------------------                    
      DO 2 I=1,NDATA                                                            
 2    WRITE(IDISK,*) ARRAY(I), QUANT(I)                                         
      STOP                                                                      
      END                                                                       
C                                                                               
      SUBROUTINE CORREL(V1,V2,NDATA,COEFF)                                      
C------------------------------------------------------------                   
C  COMPUTES CORRELATION COEFFICIENT OF TWO VECTORS V1 AND V2                    
C------------------------------------------------------------                   
      IMPLICIT REAL*8(A-H,O-Z)                                                  
      DIMENSION V1(NDATA),V2(NDATA)                                             
      SUMNUM = 0.0D0                                                            
      SUMD1  = 0.0D0                                                            
      SUMD2  = 0.0D0                                                            
      CALL XMEAN(V1,NDATA,XBAR1)                                                
      CALL XMEAN(V2,NDATA,XBAR2)                                                
      DO 1 I=1,NDATA                                                            
        SUMNUM = SUMNUM + (V1(I) - XBAR1)*(V2(I) - XBAR2)                       
        SUMD1  = SUMD1  + (V1(I) - XBAR1)**2                                    
        SUMD2  = SUMD2  + (V2(I) - XBAR2)**2                                    
 1    CONTINUE                                                                  
      COEFF = SUMNUM/(DSQRT(SUMD1*SUMD2))                                       
      RETURN                                                                    
      END                                                                       
C                                                                               
      SUBROUTINE XMEAN(ARRAY,NDATA,XBAR)                                        
      IMPLICIT REAL*8(A-H,O-Z)                                                  
      DIMENSION ARRAY(NDATA)                                                    
      XBAR = 0.D0                                                               
      DO 1 I=1,NDATA                                                            
 1    XBAR=XBAR+ARRAY(I)                                                        
      XBAR=XBAR/DFLOAT(NDATA)                                                   
      RETURN                                                                    
      END                                                                       
C                                                                               
      SUBROUTINE QQTEST(ARRAY,NDATA,IPRINT,IDISK,QUANT,PROB,IERR,RQ)            
C--------------------------------------------------------------------           
C     COMPUTES Q-Q PLOTS                                                        
C     EXTERNALS: MDNRIS IS A SUBROUTINE FROM IMSL                               
C--------------------------------------------------------------------           
      IMPLICIT REAL*8(A-H,O-Z)                                                  
      REAL*4 QNTL                                                               
      DIMENSION ARRAY(NDATA),  PROB(NDATA),  QUANT(NDATA)                       
C----------------------------------------------------                           
C     SORT DATA                                                                 
C-----------------------------------------------------                          
      CALL SORT8(ARRAY,NDATA)                                                   
C-----------------------------------------------------                          
C     COMPUTE CORRESPONDING PROBABILITIES                                       
C-----------------------------------------------------                          
      DO 1 I=1,NDATA                                                            
 1    PROB(I) = (I - 0.5D0)/ NDATA                                              
C                                                                               
C----------------------------------------------------------                     
C   CALCULATE QUANTILES OF GAUSSIAN DIASTRIBUTION                               
C---------------------------------------------------------                      
      DO 2 I=1,NDATA                                                            
      CALL MDNRIS(SNGL(PROB(I)),QNTL,IER)                                       
 2    QUANT(I) = DBLE(QNTL)                                                     
C                                                                               
C---------------------------------------------------------                      
C     CALCULATE CORRELATION COEFFICIENT                                         
C---------------------------------------------------------                      
      CALL CORREL(ARRAY,QUANT,NDATA,RQ)                                         
C---------------------------------------------------------                      
C     PRINT COMPUTED DATA                                                       
C---------------------------------------------------------                      
      WRITE(IPRINT,200)                                                         
      WRITE(IPRINT,100)                                                         
      WRITE(IPRINT,200)                                                         
      DO 4 I=1,NDATA                                                            
 4    WRITE(IPRINT,300) I, QUANT(I), ARRAY(I)                                   
      WRITE(IPRINT,200)                                                         
      WRITE(IPRINT,400) RQ                                                      
      WRITE(IPRINT,200)                                                         
 100  FORMAT(3X,'NUMBER',2X,'STANDARD NORMAL QUANTILE'                          
     #,3X,'ORDERED DATA')                                                       
 200  FORMAT('-----------------------------------------------------------       
     #----------------',/)                                                      
 300  FORMAT(' ',2X,I6,F15.4,8X,F15.4)                                          
 400  FORMAT(' ',9X,'CORRELATION COEFFICIENT=',F8.5)                            
      RETURN                                                                    
      END                                                                       
C                                                                               
      SUBROUTINE SORT8(V1,N)                                                    
      IMPLICIT REAL*8 (A-H,O-Z)                                                 
      DIMENSION V1(N)                                                           
      I = 1                                                                     
 1    I = I + I                                                                 
      IF(I.LE.N) GO TO 1                                                        
      M = I - 1                                                                 
 2    M = M / 2                                                                 
      IF(M.EQ.0) RETURN                                                         
      K = N - M                                                                 
      DO 4 J=1,K                                                                
        L = J                                                                   
 5      IF(L .LT. 1)  GO TO 4                                                   
        IF(V1(L+M) .GE. V1(L)) GO TO 4                                          
        X = V1(L+M)                                                             
        V1(L+M) = V1(L)                                                         
        V1(L)   = X                                                             
        L = L - M                                                               
        GO TO 5                                                                 
 4    CONTINUE                                                                  
      GO TO 2                                                                   
      END                                                                       
//*                                                                             
//LKED.USERLIB DD DSN=A.M12129.SELIBOBJ,DISP=SHR                                
//GO.FT10F001 DD DSN=&&TEMPPVA,DISP=(NEW,PASS),                                 
//        UNIT=DASD,SPACE=(TRK,(1,1),RLSE),                                     
//        DCB=(RECFM=FB,LRECL=80,BLKSIZE=8000)                                  
//GO.SYSIN    DD  *                                                             
    10                                                                          
      2.30                                                                      
      0.80                                                                      
      0.62                                                                      
      0.41                                                                      
      0.16                                                                      
     -0.10                                                                      
     -1.00                                                                      
      1.26                                                                      
      1.54                                                                      
      1.71                                                                      
//*                                                                             
//STEP2   EXEC  SAS,REGION=4096K                                                
//PVAFILE DD DSN=&&TEMPPVA,DISP=(OLD,DELETE)                                    
//SYSIN   DD  *                                                                 
GOPTIONS  DEVICE=ZETA HSIZE=11.0 VSIZE=8.0 BORDER;                              
*                                                                               
*------------------------------------------------------------                   
*                  INPUT DATA                                                   
*------------------------------------------------------------                   
*;                                                                              
DATA PVA;                                                                       
   INFILE PVAFILE;                                                              
   INPUT ARRAY QUANT;                                                           
*-----------------------------------------------------------;                   
*         PLOT OF QUANTILE-QUANTILE;                                            
*-----------------------------------------------------------;                   
PROC GPLOT DATA=PVA;                                                            
  PLOT ARRAY*QUANT/ VREF=0.0 HREF=0.0;                                          
  LABEL ARRAY=OBSERVED QUANTILE                                                 
        QUANT=GAUSSIAN QUANTILE;                                                
  SYMBOL V=STAR C=BLACK I=NONE;                                                 
 TITLE .C=BLACK .F=TITALIC Q-Q PLOTS FOR DEPARTURE FROM GAUSSIAN;               
//