Implement wrap conditions for simulation boxes

[seshnadathur]

In the case that the galaxies are in a simulation box with periodic boundary conditions, the standalone code should ideally respect that by implementing a wrap when writing the output positions to file.

If I haven't overlooked some edge use cases I think this could be done automatically whenever selection_function is uniform and boxsize is provided, by adding a line like positions_rec['data'] %= config['mesh']['boxsize'] in https://github.com/cosmodesi/pyrecon/blob/main/bin/pyrecon

What do you think?

[seshnadathur]

Actually there is also a use case where the data is in a simulation box so the selection function is uniform but one nevertheless still wants to provide randoms in order to obtain the shifted random positions. So perhaps instead of checking on selection_function and boxsize it would be better to just provide a option like wrap or similar for the user to set in the config as desired.

[adematti]

Yes, I agree wrap looks like a good option. Maybe one could add this option to assign_data and assign_randoms too (to wrap input positions before assigning particles to mesh). We would save this as an attribute to be used in the check here https://github.com/cosmodesi/pyrecon/blob/8be89714468d3e448089d66e09c5853809af26e9/pyrecon/iterative_fft_particle.py#L146 (which implicitly assumed that using a global line-of-sight means periodic wrapping is ok). Would you like to submit a PR for this, if you have some time?

[seshnadathur]

OK, I'll work on it later this week.

[epaillas]

Hi Sesh,

Since you are looking into this, I wanted to comment that I've run into some issues that appear to be edge cases of the periodic wrapping you were discussing. For some of my periodic box realizations, the CIC assignment fails after a few iterations. I think it could have to do with a galaxy being placed exactly at the border of the box:

[000000.00]  12-07 23:42  Main                         INFO     Loading config file ifftp.yaml.
[000000.01]  12-07 23:42  Main                         INFO     Convention is rsd.
[000000.01]  12-07 23:42  Main                         INFO     Using uniform selection function.
[000000.01]  12-07 23:42  Main                         INFO     Loading data catalog /cosma/home/analyse/epaillas/data/baojiu_recon/local_zspace/R106_S15_L2.fits.
[000000.23]  12-07 23:42  Main                         INFO     Size of catalog data is 1754711.
[000000.35]  12-07 23:42  IterativeFFTParticleReconstruction INFO     Using mesh RealMesh(boxsize=[2000. 2000. 2000.], boxcenter=[1000. 1000. 1000.], nmesh=[500 500 500], dtype=float32).
[000008.25]  12-07 23:42  IterativeFFTParticleReconstruction INFO     Running iteration 0.
[000010.96]  12-07 23:42  IterativeFFTParticleReconstruction INFO     Computing displacement field.
[000019.12]  12-07 23:43  IterativeFFTParticleReconstruction INFO     Running iteration 1.
[000027.22]  12-07 23:43  IterativeFFTParticleReconstruction INFO     Computing displacement field.
[000034.71]  12-07 23:43  IterativeFFTParticleReconstruction INFO     Running iteration 2.
[000042.39]  12-07 23:43  IterativeFFTParticleReconstruction INFO     Computing displacement field.
Index out of range: (ix,iy,iz) = (86,51,500) for (86.945,51.962,500.000)
[000064.76]  12-07 23:43  Exception                    CRITICAL 
====================================================================================================
Traceback (most recent call last):
  File "/cosma/home/analyse/epaillas/data/code/anaconda3/bin/pyrecon", line 7, in <module>
    exec(compile(f.read(), __file__, 'exec'))
  File "/cosma/home/analyse/epaillas/code/pyrecon/bin/pyrecon", line 415, in <module>
    main()
  File "/cosma/home/analyse/epaillas/code/pyrecon/bin/pyrecon", line 340, in main
    positions_rec['data'] = positions['data'] - recon.read_shifts('data', field=field)
  File "/cosma/home/analyse/epaillas/code/pyrecon/pyrecon/iterative_fft_particle.py", line 191, in read_shifts
    shifts = read_cic(self._positions_rec_data)
  File "/cosma/home/analyse/epaillas/code/pyrecon/pyrecon/iterative_fft_particle.py", line 187, in read_cic
    shifts[:,iaxis] = psi.read_cic(positions)
  File "/cosma/home/analyse/epaillas/code/pyrecon/pyrecon/mesh.py", line 561, in read_cic
    raise MeshError('Issue with read_cic')
pyrecon.mesh.MeshError: Issue with read_cic

Have you run into this problem before?

[seshnadathur]

The linked pull request should achieve this and also solve the CIC failure @epaillas mentioned.

[adematti]

Solved in PR #10. May I close this issue?

[seshnadathur]

Yes, good to close