damonge / CoLoRe

CoLoRe - Cosmological Lofty Realizations
GNU General Public License v3.0
17 stars 13 forks source link

Bad structure for 2LPT fields #63

Open cramirezpe opened 3 years ago

cramirezpe commented 3 years ago

We have been working with 2LPT CoLoRe boxes while trying to improve current DESI Ly--alpha mocks (which currently use lognormal boxes).

We've analyzed the Ly-alpha auto-correlation, the Ly-alpha-qso cross-correlation and the qso auto-correlation so far. The Ly-alpha auto does depend on the density field and it looks similar to the results obtained with the lognormal boxes.

However, when computing the Ly-alpha-qso and the qso auto, we got weird results. Therefore we wanted to inspect the differences in the quasars obtained with and without the 2LPT option set in CoLoRe.

What we observe is a few spots in the sky with a very high density of objects. I'll give a couple of examples of this. Focusing into a single pixel, let's start with the 2LPT box. In the left panel it is shown a whole healpix pixel while in the right panel it is shown a small region with an "unusual high density". The middle panel shows the location of this region inside the healpixel. 2lpt1

Then for the same configuration but just using a lognormal transformation: lognormal1

The get a better intuition of the separation between objects in both cases, I've computed distances between objects for a subsample inside the small region: histogram1

The second example shows more clearly that there is something wrong happening as the structure is box-shaped:
2lpt2 lognormal2 histogram2

The param.cfg used for this box is the following:

global = {
  prefix_out = "out"
  output_format = "FITS"
  output_density = true
  pk_filename = "/project/projectdirs/desi/users/cramirez/lya_mock_2LPT/Inputs/CoLoRe/PlanckDR12_kmax_matterpower_z0.dat"
  z_min = 1.6
  z_max = 3.79
  seed = 0
  write_pred = false
  pred_dz = 0.1
}
field_par = {
  r_smooth = 2.0
  smooth_potential = true
  n_grid = 4096
  dens_type = 2
  lpt_buffer_fraction = 0.6
  lpt_interp_type = 1
  output_lpt = 0
}
cosmo_par = {
  omega_M = 0.3147
  omega_L = 0.6853
  omega_B = 0.04904
  h = 0.6731
  w = -1.0
  ns = 0.9655
  sigma_8 = 0.830
}
srcs1 = {
  nz_filename = "/project/projectdirs/desi/users/cramirez/lya_mock_2LPT/Inputs/CoLoRe/Nz_qso_130618_2_colore1_hZs.txt"
  bias_filename = "/project/projectdirs/desi/users/cramirez/lya_mock_2LPT/Inputs/CoLoRe/Bz_qso_G18.txt"
  include_shear = false
  store_skewers = true
  gaussian_skewers = false

maybe there is something in this configuration not convenient for 2LPT that we don't understand.

And the relevant run parameters are: Screenshot 2021-01-07 022201

I tried to reproduce this same problem in smaller boxes (n_grid < 4096) without success.

Any idea of what could be happening here, or any analysis that we can perform to shed some light on this?

damonge commented 3 years ago

Hi @cramirezpe

ok, interesting. Reminds me of this problem: https://github.com/damonge/CoLoRe/issues/10

Basically, the default way in which we bias sources in colore is by doing 1+delta_qso = (1+delta_M)^b_qso Because it appears in the exponent, if b_qso is large and you get a large fluctuation in delta_M, the corresponding delta_qso can be really large, and one ends up with a few cells in the grid containing a large number of sources.

There are two ways to check if this is what's happening.

  1. Increase the smoothing scale, which will reduce the amplitude of fluctuations.
  2. Use one of the alternative bias models (which you can select in the Makefile). These other models tend to b*delta instead of delta^b for large delta, which makes things more stable.

Let me know if this doesn't help and we'll delve deeper.

andreufont commented 3 years ago

Hi @damonge - I would also bet in this direction. I think we should try running a few boxes with different schemes to populate quasars, and use the one that is "better behaved" (although not clear how to define this metric).

I'd try to not smooth the field too much, though, since it would also remove the power in the density skewers that we need for Lyman alpha. Let's start first with different biasing models.

damonge commented 3 years ago

Yes, exactly. I think the alternative bias models curbing the exponential dependence on b will be the answer here.

I guess the definition of "better behaved" depends on how you want to use these mocks. In all cases it probably means "whatever method can recover a realistic level of power down to the smallest possible scale"?

andreufont commented 3 years ago

Right, the metric should be something like that. I just didn't want to get to the point of doing a comparison to data, or doing a comparison to a particular theory (linear theory should not work that well either).

We just don't want something that is clearly wrong, like these spikes in the density field. I'd just try a couple of settings, and hopefully one will be clearly better even by eye.

I guess the linear bias model would need to be truncated to delta_q) >= -1, while the others do that by construction. The documentation does not mention the linear bias, but I saw in the code that this is BIAS_MODEL_3, and that indeed it already enforces positive densities of galaxies. @cramirezpe , are you using BIAS_MODEL_1? If so, let's try the other two and see what you get.

damonge commented 3 years ago

Yeah, the bias models are shown in common.h (search for "bias" there). BIAS_MODEL_2 should work better, I think. I also have a modified version of that one on a separate branch, in case you want to test it out