SUBROUTINE SOURCE(LIN,LOU,MPRINT,MWAVE,MTYPE,P,H,AMSOUR,PHSOUR)
C
      SAVE T,SF,AUX,PI
      DIMENSION T(3,3),SF(3),IPAR(4),PAR(6),P(3),H(3)
      COMMON/SOUR/ROS
C
C     RADIATION PATTERNS OF THE SOURCE
C
C     FOR MWAVE.EQ.0, READ THE INPUT DATA FOR RADIATION PATTERNS
C     FOR MWAVE.NE.0, COMPUTE THE RADIATION PATTERNS
C
c     READ THE INPUT DATA
c
      PI=4.*ATAN(1.)
      IF(MWAVE.NE.0) GO TO 1
      IF(MTYPE.EQ.0)RETURN
      READ(LIN,*)IPAR,PAR
      WRITE(LOU,100)IPAR,PAR
      GO TO 10
    1 CONTINUE
      IF(MTYPE.EQ.0)THEN
        AMSOUR=1.
        PHSOUR=0.
        RETURN
      END IF
      SQSL=0.
      DO 9 I=1,3
      SQSL=SQSL+P(I)*P(I)
    9 CONTINUE
   10 CONTINUE
C
C     COMPUTE THE RADIATION PATTERNS
C
      GO TO (11,21,31,41),MTYPE
C
C
C.....SINGLE FORCE
C
   11 CONTINUE
      IF(MWAVE.NE.0)GO TO 12
      AUX=PAR(2)/(4.*PI*ROS)
      SFAZ=PAR(1)
      SFDC=PAR(3)
      SF(1)=COS(SFAZ)*COS(SFDC)
      SF(2)=SIN(SFAZ)*COS(SFDC)
      SF(3)=SIN(SFDC)
      RETURN
C
   12 CONTINUE
      AMSOUR=0.
      DO 13 I=1,3
      AMSOUR=AMSOUR+H(I)*SF(I)
   13 CONTINUE
      AMSOUR=AUX*SQSL*AMSOUR
      GO TO 51
C
C.....DOUBLE COUPLE
C
   21 CONTINUE
      IF(MWAVE.NE.0)GO TO 22
      SM=PAR(2)
      FS=PAR(3)
      D=PAR(1)
      AL=PAR(4)
      SNF=SIN(FS)
      CSF=COS(FS)
      SND=SIN(D)
      CSD=COS(D)
      SNL=SIN(AL)
      CSL=COS(AL)
      A=CSL*CSF+CSD*SNL*SNF
      B=CSL*SNF-CSD*SNL*CSF
      C=SNF*SND
      D=SNL*SND
      T(1,1)=-2.*A*C
      T(1,2)=SND*(A*CSF-B*SNF)
      T(2,1)=T(1,2)
      T(1,3)=-(A*CSD-C*D)
      T(3,1)=T(1,3)
      T(2,2)=2.*SND*CSF*B
      T(2,3)=-(B*CSD+D*SND*CSF)
      T(3,2)=T(2,3)
      T(3,3)=2.*D*CSD
      IF(MPRINT.GE.1)WRITE(lou,102)T
      AUX=SM/(4.*PI*ROS)
      RETURN
C
   22 CONTINUE
      AMSOUR=0.
      DO 23 I=1,3
      DO 23 J=1,3
      AMSOUR=AMSOUR+H(I)*T(I,J)*P(J)
   23 CONTINUE
      AMSOUR=AUX*SQSL*AMSOUR
      GO TO 51
C
C.....EXPLOSIVE (IMPLOSIVE) SOURCE
C
   31 CONTINUE
      IF(MWAVE.NE.0)GO TO 32
      AUX=4.*PI*ROS
      AUX=PAR(1)*FLOAT(IPAR(1))/AUX
      RETURN
   32 CONTINUE
      AMSOUR=H(1)*P(1)+H(2)*P(2)+H(3)*P(3)
      AMSOUR=AUX*SQSL*AMSOUR
      GO TO 51
C
   41 CONTINUE
      IF(MWAVE.NE.0)GO TO 22
      T(1,1)=PAR(1)
      T(1,2)=PAR(2)
      T(2,1)=T(1,2)
      T(1,3)=PAR(3)
      T(3,1)=T(1,3)
      T(2,2)=PAR(4)
      T(2,3)=PAR(5)
      T(3,2)=T(2,3)
      T(3,3)=PAR(6)
      IF(MPRINT.GE.1)WRITE(lou,102)T
      AUX=1./(4.*PI*ROS)
      RETURN
C
   51 PHSOUR=0.
   50 IF(AMSOUR)52,53,53
   52 AMSOUR=-AMSOUR
      PHSOUR=PHSOUR-PI
C
  100 FORMAT(4I5,6F10.5)
  101 FORMAT(8F10.5)
  102 FORMAT(//2X,3(3E15.5/))
C
   53 RETURN
      END