!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!08.30.2011. BY XK @SCU!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!Input j1,j2,output the CG coefficient!!!!!!!!!!!!!!!!!!!!

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!MAIN PROGRAM!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
program main
implicit none
real    j1,j2,j3,m1,m2,m3,cg,jmax,jmin
integer m1_cont,m2_cont,j3_cont
write(*,*)"Please Input the values of j1:"
read(*,*)j1
write(*,*)"Please Input the values of j2:"
read(*,*)j2
!write(*,*)"Please Input the values of j3:"
!read(*,*)j3
!write(*,*)"Please Input the values of m1:"
!read(*,*)m1
!write(*,*)"Please Input the values of m2:"
!read(*,*)m2
!write(*,*)"Please Input the values of m3:"
!read(*,*)m3

jmax=j1+j2
jmin=abs(j1-j2)

   m1=-j1
do m1_cont=1,INT(2*j1+1)
   m2=-j2
do m2_cont=1,INT(2*j2+1)
   m3=m1+m2
   j3=max(jmin,abs(m3))
   
do j3_cont=1,INT(jmax-j3+1)
   call cgfactor(j1,j2,j3,m1,m2,m3,cg)
   write(*,*)":"
   write(*,"(1A,F4.1,1A,F4.1,1A,F4.1,1A,F4.1,1A,F4.1,1A,F4.1,1A)")"<",j1,",",j2,",",m1,",",m2,"|",j3,",",m3,">"
   write(*,*)"The Clebsh-Gordan Factor value:"
   write(*,*)cg
   j3=j3+1
end do

   m2=m2+1
end do

   m1=m1+1
end do

end program
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!SUBROUTINE ONE PROGRAM!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine cgfactor(j1,j2,j3,m1,m2,m3,cg)
implicit none
real    j1,j2,j3,m1,m2,m3
real    cgs
integer cgs1,cgs2,cgs3,cgs4,cgs5,cgs6,cgs7,cgs8,cgs9,cgs10
integer facs1,facs2,facs3,facs4,facs5,facs6,facs7,facs8,facs9,facs10,facs11
real    facsr1,facsr2,facsr3,facsr4,facsr5,facsr6,facsr7,facsr8,facsr9,facsr10,facsr11
integer vs1,vs2,vs3,vx1,vx2,vx3,vmax,vmin,v,vabs,vmod
real    cgxr,cgx,cgxv
integer cgx1,cgx2,cgx3,cgx4,cgx5,cgx6
integer facx1,facx2,facx3,facx4,facx5,facx6
real    cg

!if((m1+m2).ne.m3)then
!write(*,*)"the input values is BAD!!! Please rerun the program."
!else if((m1+m2).eq.m3) then
if((m1+m2).eq.m3) then
cgs1=j1+j2-j3
    if(cgs1.lt.0)then
      write(*,*)"The input value is BAD!!!"
!      exit
    else if(cgs1.eq.0)then
      facs1=1
    else if(cgs1.gt.0)then
      call factorial(cgs1,facs1)
    end if
    facsr1=REAL(facs1)
!write(*,*)"the facsr1 value:"
!write(*,*)facsr1
cgs2=j2+j3-j1
    if(cgs2.lt.0)then
      write(*,*)"The input value is BAD!!!"
!      exit
    else if(cgs2.eq.0)then
      facs2=1
    else if(cgs2.gt.0)then
      call factorial(cgs2,facs2)
    end if
    facsr2=REAL(facs2)
!write(*,*)"the facsr2 value:"
!write(*,*)facsr2
cgs3=j3+j1-j2
    if(cgs3.lt.0)then
      write(*,*)"The input value is BAD!!!"
!      exit
    else if(cgs3.eq.0)then
      facs3=1
    else if(cgs3.gt.0)then
      call factorial(cgs3,facs3)
    end if
    facsr3=REAL(facs3)
!write(*,*)"the facsr3 value:"
!write(*,*)facsr3
cgs4=j1+j2+j3+1
    if(cgs4.lt.0)then
      write(*,*)"The input value is BAD!!!"
!      exit
    else if(cgs4.eq.0)then
      facs4=1
    else if(cgs4.gt.0)then
      call factorial(cgs4,facs4)
    end if
    facsr4=REAL(facs4)
!write(*,*)"the facsr4 value:"
!write(*,*)facsr4
cgs5=j1+m1
    if(cgs5.lt.0)then
      write(*,*)"The input value is BAD!!!"
!      exit
    else if(cgs5.eq.0)then
      facs5=1
    else if(cgs5.gt.0)then
      call factorial(cgs5,facs5)
    end if
    facsr5=REAL(facs5)
!write(*,*)"the facsr5 value:"
!write(*,*)facsr5
cgs6=j1-m1
    if(cgs6.lt.0)then
      write(*,*)"The input value is BAD!!!"
!      exit
    else if(cgs6.eq.0)then
      facs6=1
    else if(cgs6.gt.0)then
      call factorial(cgs6,facs6)
    end if
    facsr6=REAL(facs6)
!write(*,*)"the facsr6 value:"
!write(*,*)facsr6
cgs7=j2+m2
    if(cgs7.lt.0)then
      write(*,*)"The input value is BAD!!!"
!      exit
    else if(cgs7.eq.0)then
      facs7=1
    else if(cgs7.gt.0)then
      call factorial(cgs7,facs7)
    end if
    facsr7=REAL(facs7)
!write(*,*)"the facsr7 value:"
!write(*,*)facsr7
cgs8=j2-m2
    if(cgs8.lt.0)then
      write(*,*)"The input value is BAD!!!"
!      exit
    else if(cgs8.eq.0)then
      facs8=1
    else if(cgs8.gt.0)then
      call factorial(cgs8,facs8)
    end if
    facsr8=REAL(facs8)
!write(*,*)"the facsr8 value:"
!write(*,*)facsr8
cgs9=j3+m3
    if(cgs9.lt.0)then
      write(*,*)"The input value is BAD!!!"
!      exit
    else if(cgs9.eq.0)then
      facs9=1
    else if(cgs9.gt.0)then
      call factorial(cgs9,facs9)
    end if
    facsr9=REAL(facs9)
!write(*,*)"the facsr9 value:"
!write(*,*)facsr9
cgs10=j3-m3
    if(cgs10.lt.0)then
      write(*,*)"The input value is BAD!!!"
!      exit
    else if(cgs10.eq.0)then
      facs10=1
    else if(cgs10.gt.0)then
      call factorial(cgs10,facs10)
    end if
    facsr10=REAL(facs10)
!write(*,*)"the facsr10 value:"
!write(*,*)facsr10
facs11=INT(2*j3+1)
facsr11=REAL(facs11)
!write(*,*)"the facsr11 value:"
!write(*,*)facsr11

!!!!!!!!!!!!!!!!!!!!!!!!!!!AFTER SQRTed the cgs value!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!cgs=sqrt(facsr1*facsr2*facsr3/facsr4*facsr5*facsr6*facsr7*facsr8*facsr9*facsr10*facsr11)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!NOT SQRTed the cgs value!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
cgs=facsr1*facsr2*facsr3/facsr4*facsr5*facsr6*facsr7*facsr8*facsr9*facsr10*facsr11
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

!write(*,*)"the cgs value:"
!write(*,*)cgs

vs1=j1+j2-j3
vs2=j1-m1
vs3=j2+m2
vx1=0
vx2=j1+m2-j3
vx3=j2-j3-m1

vmax=min(vs1,vs2,vs3)
vmin=max(vx1,vx2,vx3)

if(vmax.lt.vmin)then
  write(*,*)"The input value is BAD!!!"
!  exit
else if(vmax.eq.vmin)then
   v=vmax
   cgx1=v
     if(cgx1.lt.0)then
       write(*,*)"The input value is BAD!!!"
!       exit
     else if(cgx1.eq.0)then
       facx1=1
     else if(cgx1.gt.0)then
       call factorial(cgx1,facx1)
     end if
   cgx2=j1+j2-j3-v
     if(cgx2.lt.0)then
       write(*,*)"The input value is BAD!!!"
!       exit
     else if(cgx2.eq.0)then
       facx2=1
     else if(cgx2.gt.0)then
       call factorial(cgx2,facx2)
     end if
   cgx3=j1-m1-v
     if(cgx3.lt.0)then
       write(*,*)"The input value is BAD!!!"
!       exit
     else if(cgx3.eq.0)then
       facx3=1
     else if(cgx3.gt.0)then
       call factorial(cgx3,facx3)
     end if
   cgx4=j2+m2-v
     if(cgx4.lt.0)then
       write(*,*)"The input value is BAD!!!"
!       exit
     else if(cgx4.eq.0)then
       facx4=1
     else if(cgx4.gt.0)then
       call factorial(cgx4,facx4)
     end if
   cgx5=j3-j1-m2+v
     if(cgx5.lt.0)then
       write(*,*)"The input value is BAD!!!"
!       exit
     else if(cgx5.eq.0)then
       facx5=1
     else if(cgx5.gt.0)then
       call factorial(cgx5,facx5)
     end if
   cgx6=j3-j2+m1+v
     if(cgx6.lt.0)then
       write(*,*)"The input value is BAD!!!"
!       exit
     else if(cgx6.eq.0)then
       facx6=1
     else if(cgx6.gt.0)then
       call factorial(cgx6,facx6)
     end if
   
   vmod=MOD(v,2)
   
   if(vmod.eq.1)then
     cgx=1.0/REAL((-1)*facx1*facx2*facx3*facx4*facx5*facx6)
   else if(vmod.eq.0)then
     cgx=1.0/REAL(facx1*facx2*facx3*facx4*facx5*facx6)
   end if
   
else if(vmax.gt.vmin)then
cgx=0
do v=vmin,vmax,1
   cgx1=v
     if(cgx1.lt.0)then
       write(*,*)"The input value is BAD!!!"
!       exit
     else if(cgx1.eq.0)then
       facx1=1
     else if(cgx1.gt.0)then
       call factorial(cgx1,facx1)
     end if
   cgx2=j1+j2-j3-v
     if(cgx2.lt.0)then
       write(*,*)"The input value is BAD!!!"
!       exit
     else if(cgx2.eq.0)then
       facx2=1
     else if(cgx2.gt.0)then
       call factorial(cgx2,facx2)
     end if
   cgx3=j1-m1-v
     if(cgx3.lt.0)then
       write(*,*)"The input value is BAD!!!"
!       exit
     else if(cgx3.eq.0)then
       facx3=1
     else if(cgx3.gt.0)then
       call factorial(cgx3,facx3)
     end if
   cgx4=j2+m2-v
     if(cgx4.lt.0)then
       write(*,*)"The input value is BAD!!!"
!       exit
     else if(cgx4.eq.0)then
       facx4=1
     else if(cgx4.gt.0)then
       call factorial(cgx4,facx4)
     end if
   cgx5=j3-j1-m2+v
     if(cgx5.lt.0)then
       write(*,*)"The input value is BAD!!!"
!       exit
     else if(cgx5.eq.0)then
       facx5=1
     else if(cgx5.gt.0)then
       call factorial(cgx5,facx5)
     end if
   cgx6=j3-j2+m1+v
     if(cgx6.lt.0)then
       write(*,*)"The input value is BAD!!!"
!       exit
     else if(cgx6.eq.0)then
       facx6=1
     else if(cgx6.gt.0)then
       call factorial(cgx6,facx6)
     end if

   vabs=ABS(v)
   vmod=MOD(vabs,2)
   
   if(vmod.eq.1)then
     cgxv=1.0/REAL((-1)*facx1*facx2*facx3*facx4*facx5*facx6)
   else if(vmod.eq.0)then
     cgxv=1.0/REAL(facx1*facx2*facx3*facx4*facx5*facx6)
   end if
   
   cgx=cgx+cgxv
   
end do
end if
!write(*,*)"the cgx value:"
!write(*,*)cgx

!!!!!!!!!!!!!!!!!!!!!!!!!!!AFTER SQRTed the cgs value!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!if(cgx.eq.0)then
!  cg=0
!!  write(*,*)"The Clebsh-Gordan Factor value:"
!!  write(*,*)cg
!else if(cgx.ne.0)then
!  cgxr=REAL(cgx)

!!!write(*,*)"the cgxr value:"
!!!write(*,*)cgxr

!  cg=cgs*cgxr
!!  write(*,*)"The Clebsh-Gordan Factor value:"
!!  write(*,*)cg
!end if
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!NOT SQRTed the cgs value!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
if(cgx.eq.0)then
   cg=0
!  write(*,*)"The Clebsh-Gordan Factor value:"
!  write(*,*)cg
else if(cgx.gt.0)then
  cgxr=(REAL(cgx))**(2)

!!write(*,*)"the cgxr value:"
!!write(*,*)cgxr

  cg=cgs*cgxr
!  write(*,*)"The Clebsh-Gordan Factor value:"
!  write(*,*)cg
else if(cgx.lt.0)then
   cgxr=(REAL((-1)*cgx))**(2)

!!write(*,*)"the cgxr value:"
!!write(*,*)cgxr

  cg=(-1)*cgs*cgxr
!  write(*,*)"The Clebsh-Gordan Factor value:"
!  write(*,*)cg
end if
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

end if
end subroutine
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!SUBROUTINE TWO PROGRAM!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine factorial(k,fac)
implicit none
integer k,fac,fac_cont
!     if(k.lt.0)then
!     write(*,*)"the values of i, j are BAD"
!else if(k.eq.0)then
!     write(*,*)"the value of factorial:"
!     write(*,*)1
!else if(k.gt.0)then
     fac=1
     if(k.gt.0)then
     do fac_cont=1,k
        fac=fac*fac_cont
     end do
!     write(*,*)"the value of factorial:"
!     write(*,*)fac
     end if
end subroutine
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!