	FUNCTION RL(N1,L1,N2,L2,L)
	IMPLICIT REAL*8(A-H,O-Z)
C++
C TITLE: RL
C FACILITY: Subprogram of RADME, TRAMP
C AUTHOR: M.N. Harakeh, KVI, Groningen, NL.
C COMMENT:
C   In this function, multiply RL by (A**(XL/6.))
C   to get the integral in units of (FM)**L
C   A = atomic mass of the nucleus
C--
C  
	N11=N1-1
	N22=N2-1
C 
	IF (PAR(L1+L2+L).NE.1.) THEN
	   RL = 0.D0
	   RETURN
	ENDIF
	S=0.D0
C  
	DO 10 K1=1,N1
	  K11=K1-1
	DO 10 K2=1,N2
	  K22=K2-1
	  T=BINOM(N11,K11)*BINOM(N22,K22)/DFAC(2*L1+2*K11+1)
	  T=T*DFAC(2*K11+2*K22+L1+L2+L+1)*PAR(K11+K22)
	  S=S+T/DFAC(2*L2+2*K22+1)
   10	CONTINUE
C  
	S=S*SQRT(DFAC(2*L1+2*N11+1)*DFAC(2*L2+2*N22+1))
	S=S/SQRT(FAC(N11)*FAC(N22))
	RL=S/SQRT(2.D0**(N11+N22+L))
	RETURN
	END
C  
	FUNCTION ANINEJ(A1,A2,A3,B1,B2,B3,C1,C2,C3)
	IMPLICIT REAL*8(A-H,O-Z)
C++
C TITLE: ANINEJ
C FACILITY: ALLJ
C PROCEDURES CALLED: SIXJ, PAR
C AUTHOR: M.HARAKEH, KVI, GRONINGEN, NL.
C--
C  
	ANINEJ = 0.D0
C  
	X1 = MAX(ABS(A1-C3),ABS(B1-C2),ABS(A2-B3))
	X2 = MIN(A1+C3,B1+C2,A2+B3)
	N1=INT(X1)+1
	N2=INT(X2)+1
	IF (N2.LT.N1) RETURN
	AI=X1-1.
C  
	DO 10 I = N1,N2
	AI=AI+1.
	T = SIXJ(A1,B1,C1,C2,C3,AI)*SIXJ(A2,B2,C2,B1,AI,B3)
	1   *SIXJ(A3,B3,C3,AI,A1,A2)
	T = PAR(NINT(2.*AI))*(2.*AI+1.)*T
10	ANINEJ = ANINEJ+T
	RETURN
	END
C  
	FUNCTION SIXJ(A1,A2,A3,B1,B2,B3)
	IMPLICIT REAL*8(A-H,O-Z)
C++
C TITLE: SIXJ
C FACILITY: ALLJ, TRAMP
C PROCEDURES CALLED: TD
C AUTHOR: M.HARAKEH, KVI, GRONINGEN, NL.
C--
	SIXJ=PAR(NINT(A1+A2+B1+B2))*TD(A1,A2,A3,B1,B2,B3)
	RETURN
	END
C  
	FUNCTION W(A1,A2,B2,B1,A3,B3)
	IMPLICIT REAL*8(A-H,O-Z)
C++
C TITLE: W (RACAH COEFFICIENT)
C FACILITY: ALLJ
C PROCEDURES CALLED: TD
C AUTHOR: M.HARAKEH, KVI, GRONINGEN, NL.
C--
  
	W=TD(A1,A2,A3,B1,B2,B3)
	RETURN
	END
  
	FUNCTION TD(A1,A2,A3,B1,B2,B3)
	IMPLICIT REAL*8(A-H,O-Z)
C++
C TITLE: TD
C FACILITY: W, SIXJ
C PROCEDURES CALLED: BINOM, TRIREL, PAR
C AUTHOR: M.HARAKEH, KVI, GRONINGEN, NL.
C--
C  
	D(J,L,M)=1./SQRT(BINOM(J,(J+L-M)/2)*BINOM((J+L+M)/2,J)
	1*((J+L+M)/2+1))
C  
	TD = 0.D0
	IF (MOD(A1+A2+A3,1.D0).GT.1.0E-6.OR.MOD(B1+B2+A3,1.D0)
	1 .GT.1.0E-6 .OR.MOD(A1+B2+B3,1.D0).GT.1.0E-6) RETURN
C  
	J1=NINT(2.*A1)
	J2=NINT(2.*A2)
	J3=NINT(2.*A3)
	L1=NINT(2.*B1)
	L2=NINT(2.*B2)
	L3=NINT(2.*B3)
	IF (TRIREL(J1,J2,J3).LT.0..OR.TRIREL(L1,L2,J3).LT.0.
	1 .OR.TRIREL(J1,L2,L3).LT.0..OR.TRIREL(L1,J2,L3).LT.0.) RETURN
C  
	N1=IABS(MIN(-J1-L1+J3+L3,-J2-L2+J3+L3,0))
	N2=MIN(J1+J2+L1+L2+2,J1+J2-J3,L1+L2-J3,J1+L2-L3,L1+J2-L3)
	IF (N2.LT.N1) RETURN
C  
	N1=N1/2+1
	N2=N2/2+1
	SUM = 0.D0
	JL12=(J1+L1+J2+L2)/2
C  
	DO 10 I=N1,N2
	I1=I-1
	T=BINOM((J1+J2-J3)/2,I1)*BINOM((L1+L2-J3)/2,I1)/BINOM
	1(JL12+2*I1+1,3*I1)*BINOM((J1+L2-L3)/2,I1)*BINOM((L1+J2-L3)/2
	2,I1)/BINOM(3*I1,I1)*BINOM((J2-L1+L3)/2+I1,(J1+J2-J3)/2)*BINOM
	3((L1-J2+L3)/2+I1,(L1+L2-J3)/2)/BINOM(2*I1,I1)*BINOM(L3+2*I1
	4,(J2-L1+L3)/2+I1)*BINOM(JL12-L3,(J1+L2-L3)/2)*BINOM(JL12+2*I1
	5,L3+2*I1)*(JL12+2*I1+1)
   10	SUM = SUM + PAR(I1) * T
C  
	TD=SUM*D(J1,J2,J3)*D(J1,L2,L3)*D(L1,J2,L3)*D(L1,L2,J3)
	RETURN
	END
C  
	FUNCTION TJ(A1,A2,A3,B1,B2,B3)
	IMPLICIT REAL*8(A-H,O-Z)
C++
C TITLE: TJ (3J COEFFICIENT)
C FACILITY: ALLJ, TRAMP
C PROCEDURES CALLED: TS
C AUTHOR: M.HARAKEH, KVI, GRONINGEN, NL.
C--
C  
	TJ=PAR(NINT(A1-A2-B3))*TS(A1,B1,A2,B2,A3,-B3)
	RETURN
	END
C  
	FUNCTION CG(A1,B1,A2,B2,A3,B3)
	IMPLICIT REAL*8(A-H,O-Z)
C++
C TITLE: CG (CLEBSCH-GORDON COEFFICIENT)
C FACILITY: ALLJ
C PROCEDURES CALLED: TS
C AUTHOR: M.HARAKEH, KVI, GRONINGEN, NL.
C--
C  
	CG=SQRT(2.*A3+1)*TS(A1,B1,A2,B2,A3,B3)
	RETURN
	END
C  
	FUNCTION TS(A1,B1,A2,B2,A3,B3)
	IMPLICIT REAL*8(A-H,O-Z)
C++
C TITLE: TS
C FACILITY: CG, TJ
C PROCEDURES CALLED: TRIREL, BINOM, PAR
C AUTHOR: M.HARAKEH, KVI, GRONINGEN, NL.
C--
C  
	TS=0.D0
	IF (MOD(A1+A2+A3,1.D0).GT.1.0E-6) RETURN
	IF (MOD(A1-B1,1.D0).GT.1.0E-6.OR.MOD(A2-B2,1.D0).GT.
	1   1.0E-6.OR.MOD(A3-B3,1.D0).GT.1.0E-6) RETURN
	IF (MOD(2.*B3,1.D0).GT.1.0E-6) RETURN
C  
	J1=NINT(2.*A1)
	J2=NINT(2.*A2)
	J3=NINT(2.*A3)
	M1=NINT(2.*B1)
	M2=NINT(2.*B2)
	M3=NINT(2.*B3)
	IF (TRIREL(J1,J2,J3) .LT.0.) RETURN
	IF (M1+M2.NE.M3) RETURN
C  
	JM1=(J1-M1)/2
	JM2=(J2+M2)/2
	N1=IABS(MIN(J3-J2+M1,J3-J1-M2,0))
	N2=MIN(J1+J2-J3,J1-M1,J2+M2)
	IF (N2.LT.N1) RETURN
	N1=N1/2+1
	N2=N2/2+1
	S=0.D0
C  
	DO 10 I=N1,N2
	I1=I-1
	T=BINOM((J3+J1-J2)/2+I1,JM1)*BINOM((J3+J2-J1)/2+I1,JM2)/BINOM
	1(2*I1,I1)*BINOM(JM1,I1)*BINOM(JM2,I1)/BINOM(J3+2*I1,2*I1)*
	2BINOM((J1+J2-J3)/2,I1)*BINOM(J3+2*I1,(J3+J1-J2)/2+I1)
10	S=S+PAR(I1)*T
C  
	J123=(J1+J2+J3)/2
	T=BINOM(J123,J3)/(BINOM(J123,J1)*BINOM(J123,J2)*(J123+1))
	T=SQRT(T)/SQRT(BINOM(J1,(J1+M1)/2)*BINOM(J2,JM2)*BINOM
	1(J3,(J3+M3)/2))
	TS=T*S
	RETURN
	END
C  
	FUNCTION BINOM(N,M)
	IMPLICIT REAL*8(A-H,O-Z)
C++
C TITLE: BINOM
C FACILITY: ALLJ, RADME, TRAMP
C AUTHOR: M.HARAKEH, KVI, GRONINGEN, NL.
C--
C  
	DIMENSION BN(50,25)
C  
	BINOM = 1.D0
	IF (N.EQ.0 .OR. N.EQ.M .OR. M.EQ.0) RETURN
C  
	IF (N-M.GT.M) THEN
	  I = M
	ELSE
	  I = N - M
	ENDIF
C  
	IF (N.GT.50) THEN
	  J=N-I
	  AI=0.D0
	  DO 10 K=1,I
	  AI=AI+1.D0
10	  BINOM=BINOM*(J+K)/AI
	  RETURN
	ENDIF
C  
	IF (BN(N,I).EQ.0.) THEN
	  J=N-I
	  AI=0.D0
	  DO 20 K=1,I
	  AI=AI+1.D0
20	  BINOM=BINOM*(J+K)/AI
	  BN(N,I)=BINOM
	ELSE
	  BINOM = BN(N,I)
	ENDIF
	RETURN
	END
C  
	FUNCTION FAC(N)
	IMPLICIT REAL*8(A-H,O-Z)
C++
C TITLE: FAC
C FACILITY: RADME, TRAMP
C AUTHOR: M.HARAKEH, KVI, GRONINGEN, NL.
C--
C	
	DIMENSION F(50)
C  
	IF (N) 10,20,30
C  
10	FAC=0.D0
	RETURN
C  
20	FAC=1.D0
	RETURN
C  
30	IF (N.GT.50) THEN
	  FAC=1.D0
	  DO 31 I=1,N
31	  FAC=I*FAC
	  RETURN
	ENDIF
C  
	IF (F(N).EQ.0.) THEN
	  FAC=1.D0
	  DO 32 I=1,N
32	  FAC=I*FAC
	  F(N)=FAC
	ELSE
	  FAC=F(N)
	ENDIF
	RETURN
	END
C  
	FUNCTION DFAC(N)
	IMPLICIT REAL*8(A-H,O-Z)
C++
C TITLE: DFAC
C FACILITY: RADME, TRAMP
C AUTHOR: M.HARAKEH, KVI, GRONINGEN, NL.
C--
C  
	DIMENSION F(50)
C  
	IF (N) 10,20,30
C  
10 	DFAC=0.D0
	RETURN 
C  
20	DFAC=1.D0
	RETURN
C  
30	IF (PAR(N).EQ.1.) GO TO 10
	IF (N.GT.50) THEN
	  DFAC=1.D0
	  DO 31 I=1,N,2
31	  DFAC=I*DFAC
	  RETURN
	ENDIF
C  
	IF (F(N).EQ.0.) THEN
	  DFAC=1.D0
	  DO 32 I=1,N,2
32	  DFAC=I*DFAC
	  F(N)=DFAC
	ELSE
	  DFAC=F(N)
	ENDIF
	RETURN 
	END
C  
	FUNCTION TRIREL(J1,J2,J3)
	IMPLICIT REAL*8(A-H,O-Z)
C++
C TITLE: TRIREL
C FACILITY: ALLJ, TRAMP
C AUTHOR: M.HARAKEH, KVI, GRONINGEN, NL.
C--
C  
	IF (    J1+J2-J3 .LT. 0
	1 .OR.  J1-J2+J3 .LT. 0
	1 .OR. -J1+J2+J3 .LT. 0 ) THEN
	   TRIREL=-1.D0
	ELSE
	   TRIREL= 1.D0
	ENDIF
	RETURN
	END
C  
	FUNCTION PAR(N)
	IMPLICIT REAL*8(A-H,O-Z)
C--
C TITLE: PAR
C FACILITY: RADME, ALLJ, TRAMP
C AUTHOR: M.HARAKEH, KVI, GRONINGEN, NL.
C--
C  
	I=N/2
	I=2*I
	IF (I.EQ.N) THEN
	  PAR= 1.D0
	ELSE
	  PAR=-1.D0
	ENDIF
	RETURN
	END
