Lyman α radiative transfer in galaxies simulated with adaptive mesh refinement.
This is the very first attempt at publishing MoCaLaTA, the Lyα RT code I built during my Ph.D. At the moment, it is mostly meant to be used by collaborators. Hence the lack of documentation. If you find this page and would like to know more, don't hesitate to contact me.
The version of MoCaLaTA is slightly different from my "own" version, but will be updated along the way. It's serial, i.e. it isn't parallelized. Written in Fortran 95/2003 in one, big chunk, the idea is that MoCaLaTA should be virtually a plug-and-play code, requiring no installation,
make'ing or anything beyond compiling (e.g. with gfortran).
The code is free for anyone to use, copy, modify, distribute, and so on, as desribed in the License. If your work leads to a publication, I'd like you to cite
where the code is described, tested, and applied. The implementation of dust in the code is described in Laursen, Sommer-Larsen, & Andersen (2009, ApJ, 704, 1640).
MoCaLaTA has previously been used in the following publications. If you use it, please let me know and I can put you on the list.
- Laursen et al. (2019, arXiv:180607392)
- Krogager et al. (2013, MNRAS, 433, 3091)
- Laursen et al. (2013, ApJ, 766, 124)
- Jensen et al. (2013, MNRAS, 428, 1366)
- Christensen et al. (2012, MNRAS, 427, 1973)
- Noterdaeme et al. (2012, A&A, 540, 63)
- Fynbo et al. (2011, MNRAS, 413, 2481)
- Watson & Laursen (2011, A&A, 527, 104)
- Laursen et al. (2011, ApJ, 728, 52)
- Laursen et al. (2009, ApJ, 704, 1640)
MoCaLaTA assumes an adaptively refined grid, and is thus best for AMR simulations (or regular grids). If your simulation is SPH, the physical properties of the particles must first be interpolated onto a grid. One way to do this is to use the visualization tool SPLASH which can at least interpolate to a regular grid, and possibly also an AMR grid; I'll get back to that.
The AMR format assumed by MoCaLaTA is a fully threaded grid, further described in the documentation for the IGM Lyα RT code IGMtransfer. At the moment, I provide a code to convert the output from SPLASH to a (regular) grid readable by MoCaLaTA.
MoCaLaTA stands for "Monte Carlo Lyman alpha Transfer with AMR", referring to the fact that the Lyα photons are transferred in an AMR grid using a Monte Carlo scheme.
That is, photons are emitted from random locations (with the probability of being emitted from a given location proportional to the Lyα luminosity in that location) and random directions, travel a random distance through the interstellar (and intergalactic) medium (given by the optical depth), are absorbed/re-emitted by hydrogen atoms (or dust grains) in random directions (given by the phase function), thus scattering stochastically until they either escape or are absorbed by dust.
In each of the cells of MoCaLaTA's grid, the following physical parameters are needed
- HI number density
- 3D velocity
- Lyα luminosity (these data are used for calculating initial positions of photons, but are not stored during runtime)
- Metallicity (optional)
- HII number density (optional)
The metallicity and the HII density are used for calculating the dust density, which is an optional feature of MoCaLaTA.
The photons that escape form an image and a spectrum of the galaxy. The galaxy may be "observed" simultaneously from six different directions (along the Cartesian axes).
I usually use the gfortran compiler. For the auxiliary codes, you can just replace the
ifort, if that's your compiler, but for
MoCaLaTA.f90 it's a slightly different syntax. In all cases, I call the executables
Compile auxiliary codes with
gfortran -O3 splash2jesper.f90 -o splash2jesper.x gfortran -O3 buildAMR.f90 -o buildAMR.x gfortran -O3 initialEmission.f90 -o initialEmission.x
MoCaLaTA.f90 with one of the following lines
gfortran -O3 -x f95-cpp-input -Ddust -DAMR MoCaLaTA.f90 -o MoCaLaTA.x ifort -O3 -fpp -Ddust -DAMR MoCaLaTA.f90 -o MoCaLaTA.x
Assuming SPH data interpolated onto a regular grid with SPLASH, each physical parameter is written to a file in the form
do k=1,nk do j=1,nj write(*,*) (dat(i,j,k),i=1,ni) enddo enddo
A series of test SPLASH files may be created in a directory
These SPLASH data files are converted to MoCaLaTA-readable format with
./splash2jesper.x testS2J.in # Convert test data to a binary called "Jesper.bin"
In the above,
testS2J.in is an input file containing the size of the computational box, as well as the names of the SPLASH files.
test.in is the input file containing all relevant information for the "real" MoCaLaTA run — i.e. in- and output file names, gas and photon parameters, cosmology, etc.
In order for MoCaLaTA not to spend memory on storing the luminosity data, this is used first to calculate the initial position of all photons, after which the luminosity isn't used. First the grid is built, and then the positions are calculated:
./buildAMR.x < test.in # Build the grid from Jesper.bin, creating "RawData.bin" ./initialEmission.x < test.in # Calculate emission locations, creating "CellData.bin"
(notice that in this case, you shouldn't use the
As the code is running, some output is generated along the way. This can be redirected into a log file instead, using
./MoCaLaTA.x test.in > test.log
A text file containing information about each escaping photon is generated, along with up to six "IFU" cubes representing an "observation" of the galaxy along the Cartesian axes.
Analyzing the output
A Python code (under development) is provided. Currently, it can do two things
- generate an image and a spectrum from an IFU cube, and
- shows the anisotropic escape of the photons (using a Mollweide projection)
Image and spectrum
The IFU cubes are saved in the directory
SubDir/, and are called
xm.bin (for the view along the $x$ axis in the negtive direction),
xp.bin, etc. The keywords
SubDir are set in the input file (
test.in in the example), where you can also set a keyword
outlabel, which is appended to the basename of the files, just before the extension
To "observe" the galaxy in the example above from the positive direction of the $x$ axis, do (in Python):
>>> import analyzeMoCaLaTA as am >>> am.showobs('testdir','xp')
This should produce something like the following figure
To extract only the spectrum of the pixels inside an circular aperture of radius
r kpc centered on [
y0] (also in kpc), use
In addition to the IFU cubes, MoCaLaTA also spits out a text file with some information about each photon packet
WILL WRITE MORE HERE. THE FOLLOWING FIGURE WILL BE PRODUCED