C ANGSOL    SOURCE    CHAT      05/01/12    21:22:26     5004           
C   CALCUL DE L'ANGLE SOLIDE D'UN TRIEDRE
C
C    PIERRE VERPEAUX     MAI 1979
C
C
C     VERSION EN SIMPLE PRECISION   FORTRAN PUR
C
      SUBROUTINE ANGSOL(XO,XA,XB,XC,EPS,IFLAG,IFLIG)
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8 (A-H,O-Z)
-INC CCREEL
      DIMENSION XO(3),XA(3),XB(3),XC(3),XU(3),XV(3),XW(3)
C  CALCUL DES TROIS VECTEURS COMPOSANT LE TRIEDRE
      IFLIG=0
      DO 5 I=1,3
         A=XO(I)
         XU(I)=XA(I)-A
         XV(I)=XB(I)-A
         XW(I)=XC(I)-A
 5    CONTINUE
C   CALCUL DES NORMES DES VECTEURS
      U=0.D0
      V=0.D0
      W=0.D0
      DO 6 I=1,3
         U=U+XU(I)**2
         V=V+XV(I)**2
         W=W+XW(I)**2
 6    CONTINUE
      U= SQRT(U)
      V= SQRT(V)
      W= SQRT(W)
      IF(IFLAG.EQ.1) THEN
         IF (U.LT.XPETIT) GOTO 11
         IF (V.LT.XPETIT) GOTO 11
         IF (W.LT.XPETIT) GOTO 11
      ENDIF
      IF (U.LT.XPETIT) GOTO 10
      IF (V.LT.XPETIT) GOTO 10
      IF (W.LT.XPETIT) GOTO 10
C   CALCUL DU PRODUIT MIXTE DES TROIS VECTEURS
      VOL=XU(1)*(XV(2)*XW(3)-XV(3)*XW(2))-XU(2)*(XV(1)*XW(3)-XV(3)*XW(1)
     #     )+XU(3)*(XV(1)*XW(2)-XV(2)*XW(1))
      IF(IFLAG.EQ.1) THEN
         IF (ABS(VOL).LT.0.0001D0) CALL COPLAN(XO,XA,XB,XC,IFLIG)
         IF (IFLIG.EQ.1) GOTO 12
      ENDIF
C CALCUL DES ANGLES
      A= ACOS((XV(1)*XW(1)+XV(2)*XW(2)+XV(3)*XW(3))/(V*W))
      B= ACOS((XW(1)*XU(1)+XW(2)*XU(2)+XW(3)*XU(3))/(W*U))
      C= ACOS((XU(1)*XV(1)+XU(2)*XV(2)+XU(3)*XV(3))/(U*V))
      P=(A+B+C)/2.D0
C    CALCUL DES NUMERATEURS ET DENOMINATEURS DANS LE
C    SECOND MEMBRE DE LA FORMULE DE L'HUILIER
      P1=P/2.D0
      P2=(P-A)/2.D0
      P3=(P-B)/2.D0
      P4=(P-C)/2.D0
      DNUM= SIN(P1)* SIN(P2)* SIN(P3)* SIN(P4)
      IF (DNUM.LT.0.D0) DNUM=0.D0
      DNUM= SQRT(DNUM)
      DDEN= COS(P1)* COS(P2)* COS(P3)* COS(P4)
      IF (DDEN.LT.0.D0) DDEN=0.D0
      DDEN= SQRT(DDEN)
C   CALCUL DE L'ANGLE SOLIDE
      EPS= ATAN2(DNUM,DDEN)
      EPS=4.D0*EPS
      IF (VOL.LE.0.D0) EPS=-EPS
      RETURN
 10   EPS=0.D0
      RETURN
C                        UN DES VECTEURS EST NUL
 11   IFLIG=1
 12   RETURN
      END