C
C *** THE HASKELL'S MATRIX: ARS2(ARR2) ,ALS2(ALR2), 
C                           ARR1,Z(ARS1,Z), AND ALR1,Z(ALS1,Z) ***
C
      SUBROUTINE ARAL2Z(AR,ARZ,AL,ALZ,LL,DEP,OM2,CV2,WN
     *                 ,DNS,VP,VS,THK,ML,NL)
C
      IMPLICIT REAL*8(A-H,O-Z)
      DIMENSION DNS(ML),THK(ML)
      COMPLEX*16 VP(ML),VS(ML)
     *          ,AR(4,4),SAR(4,4),SBR(4,4)
     *          ,AL(2,2),SAL(2,2),SBL(2,2)
     *          ,ARZ(4,4),SBRZ(4,4)
     *          ,ALZ(2,2),SBLZ(2,2)
     *          ,WP2,WS2,CMP,CMS,BM2,PM,QM,RBM,RAM
     *          ,CQ,SQ,STQ,SDQ,CP,SP,STP,SDP
     *          ,CQZ,SQZ,CPZ,SPZ
C
      WN2=WN*WN
C
C *** INITIAL MATRIX ***
C
      DO 100 I=1,4
        DO 110 J=1,4
          AR(I,J) =(0.0D0,0.0D0)
  110   CONTINUE
        AR(I,I) =(1.0D0,0.0D0)
  100 CONTINUE
C
      DO 150 I=1,2
        DO 160 J=1,2
          AL(I,J) =(0.0D0,0.0D0)
  160   CONTINUE
        AL(I,I) =(1.0D0,0.0D0)
  150 CONTINUE
C
C *** MATRIX FROM IL=NL-1 TO LL-1 ***
C
      DO 200 IL=NL-1,LL+1,-1
C
        DO 210 I=1,4
          DO 220 J=1,4
            SAR(I,J) =AR(I,J)
  220     CONTINUE
  210   CONTINUE
C
        DO 230 I=1,2
          DO 240 J=1,2
            SAL(I,J) =AL(I,J)
  240     CONTINUE
  230   CONTINUE
C
        BM2=VS(IL)*VS(IL)
        WP2=OM2/(VP(IL)*VP(IL))
        WS2=OM2/BM2
C
        CMP=CDSQRT(WP2-WN2)
        CMS=CDSQRT(WS2-WN2)
        IF(IMAG(CMP).GT.0.D0)CMP=-CMP
        IF(IMAG(CMS).GT.0.D0)CMS=-CMS
C
        RAM=CMP/WN
        PM=CMP*THK(IL)
        CP=CDCOS(PM)
        SP=CDSIN(PM)
        STP=RAM*SP
        SDP=SP/RAM
C
        RBM=CMS/WN
        QM=CMS*THK(IL)
        CQ=CDCOS(QM)
        SQ=CDSIN(QM)
        STQ=RBM*SQ
        SDQ=SQ/RBM
C
        CALL ARALSB(SBR,SBL,CV2,DNS(IL),BM2
     *             ,CQ,STQ,SDQ,CP,STP,SDP)
C
        DO 250 I=1,4
          DO 255 J=1,4
            AR(I,J) =(0.D0,0.D0)
            DO 260 K=1,4
              AR(I,J) =AR(I,J) +SAR(I,K) *SBR(K,J)
  260       CONTINUE
  255     CONTINUE
  250   CONTINUE
C
        DO 270 I=1,2
          DO 275 J=1,2
            AL(I,J) =(0.D0,0.D0)
            DO 280 K=1,2
              AL(I,J) =AL(I,J) +SAL(I,K) *SBL(K,J)
  280       CONTINUE
  275     CONTINUE
  270   CONTINUE
C
  200 CONTINUE
C
C *** MATRIX FOR IL=LL ***
C
      DO 300 I=1,4
        DO 310 J=1,4
          SAR(I,J) =AR(I,J)
  310   CONTINUE
  300 CONTINUE
C
      DO 350 I=1,2
        DO 360 J=1,2
          SAL(I,J) =AL(I,J)
  360   CONTINUE
  350 CONTINUE
C
      BM2=VS(LL)*VS(LL)
      WP2=OM2/(VP(LL)*VP(LL))
      WS2=OM2/BM2
C
      CMP=CDSQRT(WP2-WN2)
      CMS=CDSQRT(WS2-WN2)
      IF(IMAG(CMP).GT.0.D0)CMP=-CMP
      IF(IMAG(CMS).GT.0.D0)CMS=-CMS
C
      RAM=CMP/WN
      PM=CMP*(THK(LL)-DEP)
      CP=CDCOS(PM)
      SP=CDSIN(PM)
      STP=RAM*SP
      SDP=SP/RAM
C
      RBM=CMS/WN
      QM=CMS*(THK(LL)-DEP)
      CQ=CDCOS(QM)
      SQ=CDSIN(QM)
      STQ=RBM*SQ
      SDQ=SQ/RBM
C
      CALL ARALSB(SBR,SBL,CV2,DNS(LL),BM2
     *           ,CQ,STQ,SDQ,CP,STP,SDP)
C
      CPZ= CMP*SP
      SPZ=-CMP*CP
      STP=RAM*SPZ
      SDP=SPZ/RAM
C
      CQZ= CMS*SQ
      SQZ=-CMS*CQ
      STQ=RBM*SQZ
      SDQ=SQZ/RBM
C
      CALL ARALSB(SBRZ,SBLZ,CV2,DNS(LL),BM2
     *           ,CQZ,STQ,SDQ,CPZ,STP,SDP)
C
      DO 400 I=1,4
        DO 410 J=1,4
          AR(I,J) =(0.D0,0.D0)
          ARZ(I,J)=(0.D0,0.D0)
          DO 420 K=1,4
            AR(I,J) =AR(I,J) +SAR(I,K) *SBR(K,J)
            ARZ(I,J)=ARZ(I,J)+SAR(I,K) *SBRZ(K,J)
  420     CONTINUE
  410   CONTINUE
  400 CONTINUE
C
      DO 450 I=1,2
        DO 460 J=1,2
          AL(I,J) =(0.D0,0.D0)
          ALZ(I,J)=(0.D0,0.D0)
          DO 470 K=1,2
            AL(I,J) =AL(I,J) +SAL(I,K) *SBL(K,J)
            ALZ(I,J)=ALZ(I,J)+SAL(I,K) *SBLZ(K,J)
  470     CONTINUE
  460   CONTINUE
  450 CONTINUE
C
      RETURN
      END
