Module Drawin
   !
   USE MOD_CSTES, ONLY: me => cste_m_e, mu => cste_m_u, mh => cste_m_p, Ryd => cste_Ryd, & 
                         h => cste_h, c => cste_c, kb =>  cste_k_B   
   USE MTD_STRUCT
   USE CONFIG, ONLY: Mass
   !
   Implicit None
   !
   PRIVATE
   PUBLIC :: DRAWIN_BB_BF
   !
CONTAINS
   !
   Subroutine DRAWIN_BB_BF(K, TEMP_IN, S0, S1, S2, UPSILON_TEMP_IN,I)
      Integer,                        intent(in) :: K
      Real,                           intent(in) :: TEMP_IN
      Type(ATOM)      ,               intent(in) :: S0
      Type(LEVELS)    , dimension(:), intent(in) :: S1  !N1
      Type(INDEXLINES), dimension(:), intent(in) :: S2  !N2
      Real,                          intent(out) :: UPSILON_TEMP_IN
      Integer, optional,              intent(in) :: I
      !
      Integer          :: Ni, Nj 
      Double Precision :: En_i_cm, En_j_cm, DeltaE_cm
      Real             :: loggf, Fij, gi
      Real :: massratio, dekt
      !
      If(Present(I)) Then
         !
         ! BOUND-FREE
         !
         Ni = I
         Nj = Size(S1)
         gi = S1(Ni)%g
         Fij = 1.0
      Else
         !
         ! BOUND-BOUND
         !
         Ni = S2(K)%Ni
         Nj = S2(K)%Nj
         gi = S1(Ni)%g
         loggf = S2(K)%loggf
         Fij   = 10**( loggf ) / gi
      End If
      !     
      En_i_cm   = S1(Ni)%En_cm
      En_j_cm   = S1(Nj)%En_cm 
      DeltaE_cm = En_j_cm - En_i_cm
      !
      massratio = me / mh * (mh + Mass*mu) / (mh + me)
      dekt = h * c * DeltaE_cm*100. / (kb * TEMP_IN)
      !
      UPSILON_TEMP_IN = 4*SQRT(2.)*massratio*gi*Fij*(Ryd/100.)/(DeltaE_cm)*(1+2/dekt)/dekt
      !
      !Write(*,*)"dekt,massratio,gi*Fij,UPSILON_TEMP_IN:",dekt,massratio,gi*Fij,UPSILON_TEMP_IN
      !
   ENd Subroutine DRAWIN_BB_BF
   !
End Module Drawin