Using a lens without parameter file

[FabienNugier]

I am having some trouble with insertSubstructures and more precisely the number of planes. Indeed, there is an assertion in insertSubstructures that checks the number of planes : assert(field_planes.size() == field_Nplanes_original); and this assertion fails when I call insertSubstructures on a lens that I created without the parameter file. The source of the problem seems to be that field_Nplanes_original is not initialized and thus takes a very large and unrealistic value. I tried to initialize it to zero in the constructor of the lens but then other problems appear with field_planes...

[rbmetcalf]

I think the source plane counts as a field_plane, but I don't remember exactly. We just need to go through more carefully the constructor that uses a peremeter file and make sure the other constructor does all the things it should.

[rbmetcalf]

This Lens constructor was not intended to be used by the user. There are minimal requirements for the Lens object that need to be done by hand. The original idea was that using the parameter file you would be forced to log them.

[rbmetcalf]

See PR #113

[FabienNugier]

There was a typo as the test on lengths in TreeQuad::BuildQTreeNB was made with '!=' instead of '=='. I corrected it. I also had to use setDist() in several functions: in insertSubstructures and resetSubstructure of Lens class, in order to pass the assertion on Dist == -1. In addition to that, I used setDist() in Lens::insertMainHalo, Lens::insertMainHalos, and two functions Lens::replaceMainHalos, otherwise the functions using the parameter file where also not passing the assertion Dist == -1. I am now testing both cases (with parameter file and without) in order to see if I get the same results, but it should be fine (running at the moment).

[rbmetcalf]

Well it seems like making that change helped reveal the problem.

Can you commit these changes only into this branch so that we have a version that solves this problem only.

[FabienNugier]

Yes indeed. I had not thought that this could have been the reason actually. Thanks !

I have kept the setDist() separated from the Halo construction as I guess we don't want the halo to depend on the cosmology. Am I right ?

I had done this changes into my branch but I am going to go onto the main branch and do the same modifications. Thus we don't have to merge my branch with the main one for that modification. Is that okay this way ?

[rbmetcalf]

Yes. The general philosophy was to keep the LensHalos independent of the cosmology. This is why their constructors don't take a COSMOLOGY and calculate Dist right there.
But the light-cones are all done in angular coordinates and there was a lot of projecting onto planes while preserving angular position that was causing problems before.
It became cumbersome not to have the angular position in the LensHalo and to make that consistent with the physical position requires the distance. As it is, you should be able to construct and use a LensHalo without having access to a COSMOLOGY and their action within the ray-shooting parts should not require it. It is only when you use LensHalo::getX() that the distance must have been set. So there shouldn't be any need for assert(Dist != -1) anywhere but in LensHalo::getX().

[FabienNugier]

I see. Thanks for the explanations :). I compared the output of my project with 50 trials with and without parameter file and I find the same results almost (same distributions of mu and positions in the sky). Nevertheless there is a small difference in Dist between the two cases that I have to understand. In one case I get Dist = 999.97... and the other I get 1001.39... for a redshift of zl = 0.34 and WMAP5yr cosmology. I don't understand why and I am investigating it.

[FabienNugier]

Okay I found the source of discrepancy. It is not related to the code but just because I had taken h = 0.704 in my parameter file while WMAP5yr takes h = 0.705. The two methods thus agree (i.e. with and without parameter file).

[FabienNugier]

There is a little subtlety that I just discovered. If you define the lens from the parameter file and in this parameter file you define a main halo (even if massless, like in mine), and that you create another (zero mass) main halo in the code, and insert it into the lens, then you get a problem when defining your grid as rayshooter will at some point call the main halo from the parameter file and create a Dist == -1 alert. Thus the correct way is either not to define a main halo in the parameter file (I don't know if we are allowed not to define one), or just use replaceMainHalo in order to replace the zero mass Halo of the parameter file with the one of the code (as reset will erase the parameter file one and set Dist).

[FabienNugier]

I think we can close that issue.

[rbmetcalf]

It seems like Dist should be set when the main halo is set with the parameter file.

[FabienNugier]

Indeed, I thought this could be something to impose. I did not dare to do the changes in the parameter reading routines though.