Program MieScat

  use mie_theory

  Implicit NONE

  ! -----CONSTANT-----
  ! mathematical constant
  real(8) :: PI

  ! -----INPUT-----
  ! (real) refractive index of surrounding medium
  real(8) :: Rmed
  ! refractive index of sphere = Rre + i*Rim
  real(8) :: Rre, Rim
  ! radius and wavelength  [m]
  real(8) :: rad, wl

  ! -----OUTPUT-----
  ! extinction efficiency
  real(8) :: QEXT
  ! single scattering albedo and co-albedo
  real(8) :: SSA
  real(8) :: SSCoA
  ! asymmetry factor
  real(8) :: AF


  ! size parameter
  real(8) :: X
  ! extinction cross section  [m2]
  real(8) :: SEXT

  ! -----BEGIN-----
  ! set mathematical constant
  PI = acos(-1.0D0)

  ! set input parameters
  ! (real) refractive index of surrounding medium
  Rmed = 1.00045D0
  ! refractive index of sphere = Rre + i*Rim
  Rre = 1.427D0
  Rim = 1.19D-6
  ! radius and wavelength  (same units)
  rad = 1.0D-6
  wl  = 1.0D-6

  ! output parameters
  Write(*,*) '-----input-----'
  Write(*,*) 'refractive index of surrounding medium (real)', &
    &           Rmed
  Write(*,*) 'refractive index of particle (real)          ', &
    &           Rre
  Write(*,*) 'refractive index of particle (imaginary)     ', &
    &           Rim
  Write(*,*)
  Write(*,*) 'particle radius [m]                          ', &
    &           rad
  Write(*,*) 'wavelength [m]                               ', &
    &           wl
  Write(*,*)

  ! calculation
  X = 2.0D0*PI*rad*Rmed/wl
  Call CALLBH3(                &
    & Rmed, Rre, Rim, rad, wl, &
    & QEXT, SSA, SSCoA, AF     &
    & )
  SEXT = PI*(rad*rad)*QEXT

  ! output results
  Write(*,*) '-----output-----'
  Write(*,*) 'size parameter               ', X
  Write(*,*)
  Write(*,*) 'extinction efficiency        ', QEXT
  Write(*,*) 'extinction cross section [m2]', SEXT
  Write(*,*)
  Write(*,*) 'single scattering albedo     ', SSA
  Write(*,*) 'co-single scattering albedo  ', SSCoA
  Write(*,*) 'asymmetry factor             ', AF


  Stop
End Program MieScat
