Add code to populate skewers with extra absorption systems

[andreufont]

DLAs are one of the main contaminants of large scale Lyman alpha correlations. Technically, we not only care about DLAs but also systems with lower column densities, and we usually refer to them as HCD (high-column density systems).

We could easily add random HCD following a given column density distribution, f(N,z), where N is the HI column density. However, we would like to be able to add HCD in a way that are correlated with the density field, and with the right redshift space distortions. The fact that we have the densities and the velocities in the CoLoRe skewers is one of the selling points of these mocks, because it allows us to do precisely this.

I'm adding this issue to motivate a new branch to work on adding HCDs. I believe this could fit very well within CoLoRe itself, but if you think that it is too specific to Lyman alpha we could add this somewhere else. We should be able to do this with just a single file, with a couple hundred lines.

PS: The same code could be used to add also metal absorbers, by just specifying a different bias and a different density distribution dn/dz(z).

[damonge]

@andreufont sounds like this could be done easily in CoLoRe. I guess the question is what clustering model you'd like use for these things (e.g. lognormal, LPT etc.)

[andreufont]

@damonge , @fjaviersanchez, @jfarr03

What I would do is something similar to what I did to add DLAs in the BOSS mocks, described in Font-Ribera & Miralda-Escudé (2012). In particular, take a look at sections 2.2 (that describes the column density distribution) and 2.3 (that describes the clustering).

Basically, you can start with the Gaussian field in the skewers, g(z). Then you choose a threshold nu(z) such that if you add DLAs to the cells with g(z) > nu(z) you will end up with the desired bias b(z). For now I would just use b(z) = 2.0. If we start with the Gaussian field, it is trivial to compute the threshold that we need in order to have a certain bias. We just need to remember that the Gaussian field is at z=0, so we need to set the threshold nu(z) to get the right D(z) * b(z).

Once we have the threshold nu(z), we could decide to not add DLAs to all the cells with g(z) > nu(z). We could compute, using the column density distribution f(z,N) described in section 2.2, what is the mean number of systems that we should have at a cell of a given cell, and draw a Poisson variable to decide whether we add a system or not. But then again, this could be done in post-processing outside CoLoRe, and our job here could end by providing a list of all possible positions for the DLAs. I think I would do this latter approach, since it would also allow to use the same CoLoRe run but change the column density distribution in post-process.

[andreufont]

Let's talk about it tomorrow over lunch, since I'm not 100% convinced that this should be in CoLoRe, instead of igmhub/LyaCoLoRe. If we try to add many systems, but at the same time require a high bias, we might end up having the typical problem of having to add more than one system in a given cell. This could be avoided if we added the DLAs after adding extra small scale power, in finer pixels, and this is currently done outside of CoLoRe...

[damonge]

OK, let's do that. Not sure if nu(z) would take into account the field smoothing, but you've probably thought about that :-)

[andreufont]

As I mentioned in person, I added a notebook describing how I would do that. You can find it here: https://github.com/igmhub/LyaCoLoRe/blob/master/example_scripts/threshold_dla_bias.ipynb

We also discussed that we should probably add this in igmhub/LyaCoLoRe instead of CoLoRe, so we can probably close this issue.

[andreufont]

This is now done in LyaCoLoRe, no need to add anything here.