search

cagurto / disk_envelope_model

0 stars 3 forks source link

Test 1: Bohren & Huffman code (Draine) #2

Open cagurto opened 8 years ago

cagurto commented 8 years ago

How works:

BHMIE.F is a "Mie solutions" (infinite series) to scattering, absorption and phase function of electromagnetic waves by a homogeneous sphere.

Steps:

  1. Type "make" in the directory where bhmie.f it is, which should compile bhmie.f and make_ca_cs_g.f90 and produce a code called "makeopac"
  2. Get an .lnk file from the Jena database, e.g. http://www.astro.uni-jena.de/Laboratory/OCDB/data/silicate/amorph/**pyrmg70.lnk**
  3. Edit param.inp to put the name "pyrmg70" on the first line (if you use indeed pyrmg70.lnk). The second line should be the grain size in cm. The third line the material density in gram/cm^3. And the fourth line must, for now, be "1".
  4. Now type ./makeopac and this should produce a file dustkappa_pyrmg70.inp. This will have four columns: wavelength in micron, kappa_abs in cm^2/gram, kappa_sca in cm^2/gram and the "g" anisotropy factor.

To plot the opacity file in Python:

import numpy as np import matplotlib.pyplot as plt data=np.loadtxt('dustkappa_pyrmg70.inp',skiprows=2)

plot wavelenght v/s kappa_abs in cm^2/gram

plt.plot(data[:,0],data[:,1])

plot wavelength v/s kappa_sca in cm^2/gram

plt.plot(data[:,0],data[:,2],linestyle='--') plt.xscale('log') plt.yscale('log') plt.xlabel("$\lambda\; [\mu\mathrm{m}]$") plt.ylabel("$\kappa\; [\mathrm{cm}^2/\mathrm{g}]$") figure_bohrenhuffman_pyrmg70