Hi all,
I was trying to generate light rays using yt 3.4-dev, but found that the light rays have wrong redshift intervals. For example, my simulation box is 80Mpccm/h per side, and I expect a light ray generated in a data dump at z=2 with length =1 simulation unit (i.e. 80Mpccm/h) would have a redshift interval from z=2 to z=1.92. However, the light ray generated by yt gives spans from z=2 to z=1.97.
As I look into this, I found that the function comoving_radial_distance in cosmology.py might return the wrong unit. I think the return value should be in comoving units, instead of physical units. To see this directly, I made the following change in cosmology_splice.py:
(root) ~/yt-conda/src/yt-hg/yt/analysis_modules/cosmological_observation $diff cosmology_splice.py cosmology_splice.py0 373d372 < print target_distance, distance2, z2
And it shows:
80.0 Mpccm/h 4.65320708035e+26 cm 1.97387969592 dimensionless 80.0 Mpccm/h 1.19455177073e+26 cm 1.97343368065 dimensionless 80.0 Mpccm/h 1.21507519065e+26 cm 1.97342592996 dimensionless
Then I made the following change in cosmology.py:
(root) ~/yt-conda/src/yt-hg/yt/utilities $diff cosmology.py cosmology.py0 111,112c111,112 < return self.quan((self.hubble_distance() * < trapzint(self.inverse_expansion_factor, z_i, z_f)).value,
return (self.hubble_distance() *
trapzint(self.inverse_expansion_factor, z_i,
z_f)).in_base(self.unit_system)
Then I get: 80.0 Mpccm/h 4.65320708035e+26 cmcm 1.92163908776 dimensionless 80.0 Mpccm/h 3.62771515661e+26 cmcm 1.92124702027 dimensionless
With the change of unit from cm to cmcm, the light rays have the right span.
Even though this solves my problem, I am not sure if similar problem still exists. For example, instead of making change in comoving_radial_distance, we might need to change hubble_distance into comoving units. Hopefully someone familiar with yt unit system could check this.
Thanks, Pengfei
Hi Pengfei,
Can you open an issue so we don't lose track of this?
https://bitbucket.org/yt_analysis/yt/issues/new
If you can, include a self-contained example using one of the cosmology time series datasets from yt-project.org/data so one of us can reproduce you issue locally.
Nathan
On Monday, August 8, 2016, Pengfei Chen madcpf@gmail.com wrote:
Hi all,
I was trying to generate light rays using yt 3.4-dev, but found that the light rays have wrong redshift intervals. For example, my simulation box is 80Mpccm/h per side, and I expect a light ray generated in a data dump at z=2 with length =1 simulation unit (i.e. 80Mpccm/h) would have a redshift interval from z=2 to z=1.92. However, the light ray generated by yt gives spans from z=2 to z=1.97.
As I look into this, I found that the function comoving_radial_distance in cosmology.py might return the wrong unit. I think the return value should be in comoving units, instead of physical units. To see this directly, I made the following change in cosmology_splice.py:
(root) ~/yt-conda/src/yt-hg/yt/analysis_modules/cosmological_observation $diff cosmology_splice.py cosmology_splice.py0 373d372 < print target_distance, distance2, z2
And it shows:
80.0 Mpccm/h 4.65320708035e+26 cm 1.97387969592 dimensionless 80.0 Mpccm/h 1.19455177073e+26 cm 1.97343368065 dimensionless 80.0 Mpccm/h 1.21507519065e+26 cm 1.97342592996 dimensionless
Then I made the following change in cosmology.py:
(root) ~/yt-conda/src/yt-hg/yt/utilities $diff cosmology.py cosmology.py0 111,112c111,112 < return self.quan((self.hubble_distance() * < trapzint(self.inverse_expansion_factor, z_i,
return (self.hubble_distance() *
trapzint(self.inverse_expansion_factor, z_i,
z_f)).in_base(self.unit_system)
Then I get: 80.0 Mpccm/h 4.65320708035e+26 cmcm 1.92163908776 dimensionless 80.0 Mpccm/h 3.62771515661e+26 cmcm 1.92124702027 dimensionless
With the change of unit from cm to cmcm, the light rays have the right span.
Even though this solves my problem, I am not sure if similar problem still exists. For example, instead of making change in comoving_radial_distance, we might need to change hubble_distance into comoving units. Hopefully someone familiar with yt unit system could check this.
Thanks, Pengfei
Hi Nathan,
Thanks for your reply! Please let me interact with Britton to make sure if it's not my error. Then I'll open an issue.
Thanks, Pengfei
On Mon, Aug 8, 2016 at 7:06 AM, Nathan Goldbaum nathan12343@gmail.com wrote:
Hi Pengfei,
Can you open an issue so we don't lose track of this?
https://bitbucket.org/yt_analysis/yt/issues/new
If you can, include a self-contained example using one of the cosmology time series datasets from yt-project.org/data so one of us can reproduce you issue locally.
Nathan
On Monday, August 8, 2016, Pengfei Chen madcpf@gmail.com wrote:
Hi all,
I was trying to generate light rays using yt 3.4-dev, but found that the light rays have wrong redshift intervals. For example, my simulation box is 80Mpccm/h per side, and I expect a light ray generated in a data dump at z=2 with length =1 simulation unit (i.e. 80Mpccm/h) would have a redshift interval from z=2 to z=1.92. However, the light ray generated by yt gives spans from z=2 to z=1.97.
As I look into this, I found that the function comoving_radial_distance in cosmology.py might return the wrong unit. I think the return value should be in comoving units, instead of physical units. To see this directly, I made the following change in cosmology_splice.py:
(root) ~/yt-conda/src/yt-hg/yt/analysis_modules/cosmological_observation $diff cosmology_splice.py cosmology_splice.py0 373d372 < print target_distance, distance2, z2
And it shows:
80.0 Mpccm/h 4.65320708035e+26 cm 1.97387969592 dimensionless 80.0 Mpccm/h 1.19455177073e+26 cm 1.97343368065 dimensionless 80.0 Mpccm/h 1.21507519065e+26 cm 1.97342592996 dimensionless
Then I made the following change in cosmology.py:
(root) ~/yt-conda/src/yt-hg/yt/utilities $diff cosmology.py cosmology.py0 111,112c111,112 < return self.quan((self.hubble_distance() * < trapzint(self.inverse_expansion_factor, z_i,
return (self.hubble_distance() *
trapzint(self.inverse_expansion_factor, z_i,
z_f)).in_base(self.unit_system)
Then I get: 80.0 Mpccm/h 4.65320708035e+26 cmcm 1.92163908776 dimensionless 80.0 Mpccm/h 3.62771515661e+26 cmcm 1.92124702027 dimensionless
With the change of unit from cm to cmcm, the light rays have the right span.
Even though this solves my problem, I am not sure if similar problem still exists. For example, instead of making change in comoving_radial_distance, we might need to change hubble_distance into comoving units. Hopefully someone familiar with yt unit system could check this.
Thanks, Pengfei
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
Hi Pengfei,
I'm not able to reproduce what you're seeing, but before I go into that, I can explain a few things about how the cosmology calculator works. In yt, we are only able to define a comoving unit system in relation to a proper unit system at a specific redshift. We are unable to define a comoving unit system that is the base unit system. This is a problem for the cosmology calculator, which only calculates comoving quantities. Therefore, in this case, we are forced to set the comoving and proper unit systems to be the same. You can see this for yourself by running the following script:
from yt.utilities.cosmology import Cosmology co = Cosmology() print (co.comoving_radial_distance(1.92, 2).to("Mpc/h")) print (co.comoving_radial_distance(1.92, 2).to("Mpccm/h"))
Both will be the same answer.
Now, back to the issue you're seeing. I ran the following script: http://paste.yt-project.org/show/6741/ to calculate the redshifts dumps needed for a 80 Mpc/h box to go from z = 2 to z = 1.92 and got the following result: CosmologyOutputRedshift[0] = 2.000 CosmologyOutputRedshift[1] = 1.926
This is pretty close to what you expected. Additionally, if I run the following script, http://paste.yt-project.org/show/6743/ with the above redshifts added to the parameter file, I get that dz_max (the change in redshift when moving 80 Mpccm/h at z = 0) of 0.0746521607364, which also seems to agree with what you are finding.
Can you share with me the simulation parameter file and script that you're using to make your light rays? I can take a look at that and see if I still get what I expect.
Britton
On Mon, Aug 8, 2016 at 6:13 AM, Pengfei Chen madcpf@gmail.com wrote:
Hi all,
I was trying to generate light rays using yt 3.4-dev, but found that the light rays have wrong redshift intervals. For example, my simulation box is 80Mpccm/h per side, and I expect a light ray generated in a data dump at z=2 with length =1 simulation unit (i.e. 80Mpccm/h) would have a redshift interval from z=2 to z=1.92. However, the light ray generated by yt gives spans from z=2 to z=1.97.
As I look into this, I found that the function comoving_radial_distance in cosmology.py might return the wrong unit. I think the return value should be in comoving units, instead of physical units. To see this directly, I made the following change in cosmology_splice.py:
(root) ~/yt-conda/src/yt-hg/yt/analysis_modules/cosmological_observation $diff cosmology_splice.py cosmology_splice.py0 373d372 < print target_distance, distance2, z2
And it shows:
80.0 Mpccm/h 4.65320708035e+26 cm 1.97387969592 dimensionless 80.0 Mpccm/h 1.19455177073e+26 cm 1.97343368065 dimensionless 80.0 Mpccm/h 1.21507519065e+26 cm 1.97342592996 dimensionless
Then I made the following change in cosmology.py:
(root) ~/yt-conda/src/yt-hg/yt/utilities $diff cosmology.py cosmology.py0 111,112c111,112 < return self.quan((self.hubble_distance() * < trapzint(self.inverse_expansion_factor, z_i,
return (self.hubble_distance() *
trapzint(self.inverse_expansion_factor, z_i,
z_f)).in_base(self.unit_system)
Then I get: 80.0 Mpccm/h 4.65320708035e+26 cmcm 1.92163908776 dimensionless 80.0 Mpccm/h 3.62771515661e+26 cmcm 1.92124702027 dimensionless
With the change of unit from cm to cmcm, the light rays have the right span.
Even though this solves my problem, I am not sure if similar problem still exists. For example, instead of making change in comoving_radial_distance, we might need to change hubble_distance into comoving units. Hopefully someone familiar with yt unit system could check this.
Thanks, Pengfei
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
Hi Britton,
Thanks very much for your test! The error happens when I make light rays from one single dataset. I tried to reproduce the error with the dataset under src/yt-hg/tests. The code I used is here:
http://paste.yt-project.org/show/6751/ After running the code, I did
h5ls -d lightray.h5/grid/redshift and found the last redshift was 9.98999107782409. However, the correct value should be ~ 9.979.
When I call LightRay this way, it will call the _deltaz_forward function in cosmology_splice.py, which will call the comoving_radial_distance in cosmology.py:
distance2 = self.cosmology.comoving_radial_distance(z2, z)
, where distance2 is not in comoving units. Then it will calculate z2:
z2 = ((target_distance - distance2) / m.in_units("Mpccm / h"))
. I think the subtraction here is the problem, since it tries to convert distance2 to comoving units, even though the numerical value of distance2 is already comoving. So a wrong factor of (1+z) is multiplied.
Any comment is appreciated.
Thanks again for your help! Pengfei
On Mon, Aug 8, 2016 at 7:10 AM, Britton Smith brittonsmith@gmail.com wrote:
Hi Pengfei,
I'm not able to reproduce what you're seeing, but before I go into that, I can explain a few things about how the cosmology calculator works. In yt, we are only able to define a comoving unit system in relation to a proper unit system at a specific redshift. We are unable to define a comoving unit system that is the base unit system. This is a problem for the cosmology calculator, which only calculates comoving quantities. Therefore, in this case, we are forced to set the comoving and proper unit systems to be the same. You can see this for yourself by running the following script:
from yt.utilities.cosmology import Cosmology co = Cosmology() print (co.comoving_radial_distance(1.92, 2).to("Mpc/h")) print (co.comoving_radial_distance(1.92, 2).to("Mpccm/h"))
Both will be the same answer.
Now, back to the issue you're seeing. I ran the following script: http://paste.yt-project.org/show/6741/ to calculate the redshifts dumps needed for a 80 Mpc/h box to go from z = 2 to z = 1.92 and got the following result: CosmologyOutputRedshift[0] = 2.000 CosmologyOutputRedshift[1] = 1.926
This is pretty close to what you expected. Additionally, if I run the following script, http://paste.yt-project.org/show/6743/ with the above redshifts added to the parameter file, I get that dz_max (the change in redshift when moving 80 Mpccm/h at z = 0) of 0.0746521607364, which also seems to agree with what you are finding.
Can you share with me the simulation parameter file and script that you're using to make your light rays? I can take a look at that and see if I still get what I expect.
Britton
On Mon, Aug 8, 2016 at 6:13 AM, Pengfei Chen madcpf@gmail.com wrote:
Hi all,
I was trying to generate light rays using yt 3.4-dev, but found that the light rays have wrong redshift intervals. For example, my simulation box is 80Mpccm/h per side, and I expect a light ray generated in a data dump at z=2 with length =1 simulation unit (i.e. 80Mpccm/h) would have a redshift interval from z=2 to z=1.92. However, the light ray generated by yt gives spans from z=2 to z=1.97.
As I look into this, I found that the function comoving_radial_distance in cosmology.py might return the wrong unit. I think the return value should be in comoving units, instead of physical units. To see this directly, I made the following change in cosmology_splice.py:
(root) ~/yt-conda/src/yt-hg/yt/analysis_modules/cosmological_observation $diff cosmology_splice.py cosmology_splice.py0 373d372 < print target_distance, distance2, z2
And it shows:
80.0 Mpccm/h 4.65320708035e+26 cm 1.97387969592 dimensionless 80.0 Mpccm/h 1.19455177073e+26 cm 1.97343368065 dimensionless 80.0 Mpccm/h 1.21507519065e+26 cm 1.97342592996 dimensionless
Then I made the following change in cosmology.py:
(root) ~/yt-conda/src/yt-hg/yt/utilities $diff cosmology.py cosmology.py0 111,112c111,112 < return self.quan((self.hubble_distance() * < trapzint(self.inverse_expansion_factor, z_i,
return (self.hubble_distance() *
trapzint(self.inverse_expansion_factor, z_i,
z_f)).in_base(self.unit_system)
Then I get: 80.0 Mpccm/h 4.65320708035e+26 cmcm 1.92163908776 dimensionless 80.0 Mpccm/h 3.62771515661e+26 cmcm 1.92124702027 dimensionless
With the change of unit from cm to cmcm, the light rays have the right span.
Even though this solves my problem, I am not sure if similar problem still exists. For example, instead of making change in comoving_radial_distance, we might need to change hubble_distance into comoving units. Hopefully someone familiar with yt unit system could check this.
Thanks, Pengfei
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
Hi Pengfei,
Firstly, Nathan, JohnZ, Matt, Cameron, if you could have a look at the last paragraph of this response, I would appreciate it.
Ah, I see now. I actually get a different answer for the redshift interval when I do it by hand. I get 9.994 instead of 9.979. Still, our answers, differ from the value calculated by the code for the reason you describe, essentially that the comoving to proper conversion is different between a cosmology object created from a dataset (where a specific redshift is defined) and one created by hand. In the one created from a dataset, the comoving and proper frame differ by the factor of (1+z), whereas in the one created by hand, they are the same for the reasons I described in my previous email.
This is definitely a problem. We cannot be doing cosmology distance calculations in two different unit systems, one where comoving and proper are the same and one where they are different. As strange as this may sound, the one that is least correct is using the system where the comoving and proper frames are different. That's because something like a comoving radial distance can only be expressed in the comoving frame and has no direct translation to the proper frame via some factor of (1 + z).
The solution to this, as I see it, is the following: The cosmology object that gets associated with a dataset must have its unit_registry overridden to set the comoving and proper frames equal to each other. A better solution would be to even remove the proper frame entirely from the cosmology object's unit_registry, though I don't even know if that's possible. I would appreciate the input of someone with more intimate knowledge of the unit system, like Nathan, Matt, or John Z. As a quick test, I experimented with the case where we create a cosmology object from a dataset by passing a copy of the unit_registry and then modifying the comoving frame. I ran into some difficulties because it looks like, for example, Mpc is define independently of pc, and so I would have to modify both pccm and Mpcm, and perhaps a number of others that I am unaware of. If anyone has ideas on how best to do this, please say so.
In conclusion, I do think there is an issue that needs to be fixed. Input here would be appreciated. Thanks, Pengfei, for identifying this issue!
Britton
On Mon, Aug 8, 2016 at 10:25 PM, Pengfei Chen madcpf@gmail.com wrote:
Hi Britton,
Thanks very much for your test! The error happens when I make light rays from one single dataset. I tried to reproduce the error with the dataset under src/yt-hg/tests. The code I used is here:
http://paste.yt-project.org/show/6751/ After running the code, I did
h5ls -d lightray.h5/grid/redshift and found the last redshift was 9.98999107782409. However, the correct value should be ~ 9.979.
When I call LightRay this way, it will call the _deltaz_forward function in cosmology_splice.py, which will call the comoving_radial_distance in cosmology.py:
distance2 = self.cosmology.comoving_radial_distance(z2, z)
, where distance2 is not in comoving units. Then it will calculate z2:
z2 = ((target_distance - distance2) / m.in_units("Mpccm / h"))
. I think the subtraction here is the problem, since it tries to convert distance2 to comoving units, even though the numerical value of distance2 is already comoving. So a wrong factor of (1+z) is multiplied.
Any comment is appreciated.
Thanks again for your help! Pengfei
On Mon, Aug 8, 2016 at 7:10 AM, Britton Smith brittonsmith@gmail.com wrote:
Hi Pengfei,
I'm not able to reproduce what you're seeing, but before I go into that, I can explain a few things about how the cosmology calculator works. In yt, we are only able to define a comoving unit system in relation to a proper unit system at a specific redshift. We are unable to define a comoving unit system that is the base unit system. This is a problem for the cosmology calculator, which only calculates comoving quantities. Therefore, in this case, we are forced to set the comoving and proper unit systems to be the same. You can see this for yourself by running the following script:
from yt.utilities.cosmology import Cosmology co = Cosmology() print (co.comoving_radial_distance(1.92, 2).to("Mpc/h")) print (co.comoving_radial_distance(1.92, 2).to("Mpccm/h"))
Both will be the same answer.
Now, back to the issue you're seeing. I ran the following script: http://paste.yt-project.org/show/6741/ to calculate the redshifts dumps needed for a 80 Mpc/h box to go from z = 2 to z = 1.92 and got the following result: CosmologyOutputRedshift[0] = 2.000 CosmologyOutputRedshift[1] = 1.926
This is pretty close to what you expected. Additionally, if I run the following script, http://paste.yt-project.org/show/6743/ with the above redshifts added to the parameter file, I get that dz_max (the change in redshift when moving 80 Mpccm/h at z = 0) of 0.0746521607364, which also seems to agree with what you are finding.
Can you share with me the simulation parameter file and script that you're using to make your light rays? I can take a look at that and see if I still get what I expect.
Britton
On Mon, Aug 8, 2016 at 6:13 AM, Pengfei Chen madcpf@gmail.com wrote:
Hi all,
I was trying to generate light rays using yt 3.4-dev, but found that the light rays have wrong redshift intervals. For example, my simulation box is 80Mpccm/h per side, and I expect a light ray generated in a data dump at z=2 with length =1 simulation unit (i.e. 80Mpccm/h) would have a redshift interval from z=2 to z=1.92. However, the light ray generated by yt gives spans from z=2 to z=1.97.
As I look into this, I found that the function comoving_radial_distance in cosmology.py might return the wrong unit. I think the return value should be in comoving units, instead of physical units. To see this directly, I made the following change in cosmology_splice.py:
(root) ~/yt-conda/src/yt-hg/yt/analysis_modules/cosmological_observation $diff cosmology_splice.py cosmology_splice.py0 373d372 < print target_distance, distance2, z2
And it shows:
80.0 Mpccm/h 4.65320708035e+26 cm 1.97387969592 dimensionless 80.0 Mpccm/h 1.19455177073e+26 cm 1.97343368065 dimensionless 80.0 Mpccm/h 1.21507519065e+26 cm 1.97342592996 dimensionless
Then I made the following change in cosmology.py:
(root) ~/yt-conda/src/yt-hg/yt/utilities $diff cosmology.py cosmology.py0 111,112c111,112 < return self.quan((self.hubble_distance() * < trapzint(self.inverse_expansion_factor, z_i,
return (self.hubble_distance() *
trapzint(self.inverse_expansion_factor, z_i,
z_f)).in_base(self.unit_system)
Then I get: 80.0 Mpccm/h 4.65320708035e+26 cmcm 1.92163908776 dimensionless 80.0 Mpccm/h 3.62771515661e+26 cmcm 1.92124702027 dimensionless
With the change of unit from cm to cmcm, the light rays have the right span.
Even though this solves my problem, I am not sure if similar problem still exists. For example, instead of making change in comoving_radial_distance, we might need to change hubble_distance into comoving units. Hopefully someone familiar with yt unit system could check this.
Thanks, Pengfei
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
Can one of you file an issue that summarizes the problem so we don't lose track?
On Tuesday, August 9, 2016, Britton Smith brittonsmith@gmail.com wrote:
Hi Pengfei,
Firstly, Nathan, JohnZ, Matt, Cameron, if you could have a look at the last paragraph of this response, I would appreciate it.
Ah, I see now. I actually get a different answer for the redshift interval when I do it by hand. I get 9.994 instead of 9.979. Still, our answers, differ from the value calculated by the code for the reason you describe, essentially that the comoving to proper conversion is different between a cosmology object created from a dataset (where a specific redshift is defined) and one created by hand. In the one created from a dataset, the comoving and proper frame differ by the factor of (1+z), whereas in the one created by hand, they are the same for the reasons I described in my previous email.
This is definitely a problem. We cannot be doing cosmology distance calculations in two different unit systems, one where comoving and proper are the same and one where they are different. As strange as this may sound, the one that is least correct is using the system where the comoving and proper frames are different. That's because something like a comoving radial distance can only be expressed in the comoving frame and has no direct translation to the proper frame via some factor of (1 + z).
The solution to this, as I see it, is the following: The cosmology object that gets associated with a dataset must have its unit_registry overridden to set the comoving and proper frames equal to each other. A better solution would be to even remove the proper frame entirely from the cosmology object's unit_registry, though I don't even know if that's possible. I would appreciate the input of someone with more intimate knowledge of the unit system, like Nathan, Matt, or John Z. As a quick test, I experimented with the case where we create a cosmology object from a dataset by passing a copy of the unit_registry and then modifying the comoving frame. I ran into some difficulties because it looks like, for example, Mpc is define independently of pc, and so I would have to modify both pccm and Mpcm, and perhaps a number of others that I am unaware of. If anyone has ideas on how best to do this, please say so.
In conclusion, I do think there is an issue that needs to be fixed. Input here would be appreciated. Thanks, Pengfei, for identifying this issue!
Britton
On Mon, Aug 8, 2016 at 10:25 PM, Pengfei Chen <madcpf@gmail.com
Hi Britton,
Thanks very much for your test! The error happens when I make light rays from one single dataset. I tried to reproduce the error with the dataset under src/yt-hg/tests. The code I used is here:
http://paste.yt-project.org/show/6751/ After running the code, I did
h5ls -d lightray.h5/grid/redshift and found the last redshift was 9.98999107782409. However, the correct value should be ~ 9.979.
When I call LightRay this way, it will call the _deltaz_forward function in cosmology_splice.py, which will call the comoving_radial_distance in cosmology.py:
distance2 = self.cosmology.comoving_radial_distance(z2, z)
, where distance2 is not in comoving units. Then it will calculate z2:
z2 = ((target_distance - distance2) / m.in_units("Mpccm /
h")) + z2
. I think the subtraction here is the problem, since it tries to convert distance2 to comoving units, even though the numerical value of distance2 is already comoving. So a wrong factor of (1+z) is multiplied.
Any comment is appreciated.
Thanks again for your help! Pengfei
On Mon, Aug 8, 2016 at 7:10 AM, Britton Smith <brittonsmith@gmail.com
Hi Pengfei,
I'm not able to reproduce what you're seeing, but before I go into that, I can explain a few things about how the cosmology calculator works. In yt, we are only able to define a comoving unit system in relation to a proper unit system at a specific redshift. We are unable to define a comoving unit system that is the base unit system. This is a problem for the cosmology calculator, which only calculates comoving quantities. Therefore, in this case, we are forced to set the comoving and proper unit systems to be the same. You can see this for yourself by running the following script:
from yt.utilities.cosmology import Cosmology co = Cosmology() print (co.comoving_radial_distance(1.92, 2).to("Mpc/h")) print (co.comoving_radial_distance(1.92, 2).to("Mpccm/h"))
Both will be the same answer.
Now, back to the issue you're seeing. I ran the following script: http://paste.yt-project.org/show/6741/ to calculate the redshifts dumps needed for a 80 Mpc/h box to go from z = 2 to z = 1.92 and got the following result: CosmologyOutputRedshift[0] = 2.000 CosmologyOutputRedshift[1] = 1.926
This is pretty close to what you expected. Additionally, if I run the following script, http://paste.yt-project.org/show/6743/ with the above redshifts added to the parameter file, I get that dz_max (the change in redshift when moving 80 Mpccm/h at z = 0) of 0.0746521607364, which also seems to agree with what you are finding.
Can you share with me the simulation parameter file and script that you're using to make your light rays? I can take a look at that and see if I still get what I expect.
Britton
On Mon, Aug 8, 2016 at 6:13 AM, Pengfei Chen <madcpf@gmail.com
Hi all,
I was trying to generate light rays using yt 3.4-dev, but found that the light rays have wrong redshift intervals. For example, my simulation box is 80Mpccm/h per side, and I expect a light ray generated in a data dump at z=2 with length =1 simulation unit (i.e. 80Mpccm/h) would have a redshift interval from z=2 to z=1.92. However, the light ray generated by yt gives spans from z=2 to z=1.97.
As I look into this, I found that the function comoving_radial_distance in cosmology.py might return the wrong unit. I think the return value should be in comoving units, instead of physical units. To see this directly, I made the following change in cosmology_splice.py:
(root) ~/yt-conda/src/yt-hg/yt/analysis_modules/cosmological_observation $diff cosmology_splice.py cosmology_splice.py0 373d372 < print target_distance, distance2, z2
And it shows:
80.0 Mpccm/h 4.65320708035e+26 cm 1.97387969592 dimensionless 80.0 Mpccm/h 1.19455177073e+26 cm 1.97343368065 dimensionless 80.0 Mpccm/h 1.21507519065e+26 cm 1.97342592996 dimensionless
Then I made the following change in cosmology.py:
(root) ~/yt-conda/src/yt-hg/yt/utilities $diff cosmology.py cosmology.py0 111,112c111,112 < return self.quan((self.hubble_distance() * < trapzint(self.inverse_expansion_factor, z_i,
return (self.hubble_distance() *
trapzint(self.inverse_expansion_factor, z_i,
z_f)).in_base(self.unit_system)
Then I get: 80.0 Mpccm/h 4.65320708035e+26 cmcm 1.92163908776 dimensionless 80.0 Mpccm/h 3.62771515661e+26 cmcm 1.92124702027 dimensionless
With the change of unit from cm to cmcm, the light rays have the right span.
Even though this solves my problem, I am not sure if similar problem still exists. For example, instead of making change in comoving_radial_distance, we might need to change hubble_distance into comoving units. Hopefully someone familiar with yt unit system could check this.
Thanks, Pengfei
yt-users mailing list yt-users@lists.spacepope.org
http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
yt-users mailing list yt-users@lists.spacepope.org
http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
yt-users mailing list yt-users@lists.spacepope.org
http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
Hi Nathan,
Thanks for your reminding. I just opened an issue: https://bitbucket.org/yt_analysis/yt/issues/1258/potential-unit-error-in-cos...
Thanks, Pengfei
On Tue, Aug 9, 2016 at 5:02 AM, Nathan Goldbaum nathan12343@gmail.com wrote:
Can one of you file an issue that summarizes the problem so we don't lose track?
On Tuesday, August 9, 2016, Britton Smith brittonsmith@gmail.com wrote:
Hi Pengfei,
Firstly, Nathan, JohnZ, Matt, Cameron, if you could have a look at the last paragraph of this response, I would appreciate it.
Ah, I see now. I actually get a different answer for the redshift interval when I do it by hand. I get 9.994 instead of 9.979. Still, our answers, differ from the value calculated by the code for the reason you describe, essentially that the comoving to proper conversion is different between a cosmology object created from a dataset (where a specific redshift is defined) and one created by hand. In the one created from a dataset, the comoving and proper frame differ by the factor of (1+z), whereas in the one created by hand, they are the same for the reasons I described in my previous email.
This is definitely a problem. We cannot be doing cosmology distance calculations in two different unit systems, one where comoving and proper are the same and one where they are different. As strange as this may sound, the one that is least correct is using the system where the comoving and proper frames are different. That's because something like a comoving radial distance can only be expressed in the comoving frame and has no direct translation to the proper frame via some factor of (1 + z).
The solution to this, as I see it, is the following: The cosmology object that gets associated with a dataset must have its unit_registry overridden to set the comoving and proper frames equal to each other. A better solution would be to even remove the proper frame entirely from the cosmology object's unit_registry, though I don't even know if that's possible. I would appreciate the input of someone with more intimate knowledge of the unit system, like Nathan, Matt, or John Z. As a quick test, I experimented with the case where we create a cosmology object from a dataset by passing a copy of the unit_registry and then modifying the comoving frame. I ran into some difficulties because it looks like, for example, Mpc is define independently of pc, and so I would have to modify both pccm and Mpcm, and perhaps a number of others that I am unaware of. If anyone has ideas on how best to do this, please say so.
In conclusion, I do think there is an issue that needs to be fixed. Input here would be appreciated. Thanks, Pengfei, for identifying this issue!
Britton
On Mon, Aug 8, 2016 at 10:25 PM, Pengfei Chen madcpf@gmail.com wrote:
Hi Britton,
Thanks very much for your test! The error happens when I make light rays from one single dataset. I tried to reproduce the error with the dataset under src/yt-hg/tests. The code I used is here:
http://paste.yt-project.org/show/6751/ After running the code, I did
h5ls -d lightray.h5/grid/redshift and found the last redshift was 9.98999107782409. However, the correct value should be ~ 9.979.
When I call LightRay this way, it will call the _deltaz_forward function in cosmology_splice.py, which will call the comoving_radial_distance in cosmology.py:
distance2 = self.cosmology.comoving_radial_distance(z2, z)
, where distance2 is not in comoving units. Then it will calculate z2:
z2 = ((target_distance - distance2) / m.in_units("Mpccm /
h")) + z2
. I think the subtraction here is the problem, since it tries to convert distance2 to comoving units, even though the numerical value of distance2 is already comoving. So a wrong factor of (1+z) is multiplied.
Any comment is appreciated.
Thanks again for your help! Pengfei
On Mon, Aug 8, 2016 at 7:10 AM, Britton Smith brittonsmith@gmail.com wrote:
Hi Pengfei,
I'm not able to reproduce what you're seeing, but before I go into that, I can explain a few things about how the cosmology calculator works. In yt, we are only able to define a comoving unit system in relation to a proper unit system at a specific redshift. We are unable to define a comoving unit system that is the base unit system. This is a problem for the cosmology calculator, which only calculates comoving quantities. Therefore, in this case, we are forced to set the comoving and proper unit systems to be the same. You can see this for yourself by running the following script:
from yt.utilities.cosmology import Cosmology co = Cosmology() print (co.comoving_radial_distance(1.92, 2).to("Mpc/h")) print (co.comoving_radial_distance(1.92, 2).to("Mpccm/h"))
Both will be the same answer.
Now, back to the issue you're seeing. I ran the following script: http://paste.yt-project.org/show/6741/ to calculate the redshifts dumps needed for a 80 Mpc/h box to go from z = 2 to z = 1.92 and got the following result: CosmologyOutputRedshift[0] = 2.000 CosmologyOutputRedshift[1] = 1.926
This is pretty close to what you expected. Additionally, if I run the following script, http://paste.yt-project.org/show/6743/ with the above redshifts added to the parameter file, I get that dz_max (the change in redshift when moving 80 Mpccm/h at z = 0) of 0.0746521607364, which also seems to agree with what you are finding.
Can you share with me the simulation parameter file and script that you're using to make your light rays? I can take a look at that and see if I still get what I expect.
Britton
On Mon, Aug 8, 2016 at 6:13 AM, Pengfei Chen madcpf@gmail.com wrote:
Hi all,
I was trying to generate light rays using yt 3.4-dev, but found that the light rays have wrong redshift intervals. For example, my simulation box is 80Mpccm/h per side, and I expect a light ray generated in a data dump at z=2 with length =1 simulation unit (i.e. 80Mpccm/h) would have a redshift interval from z=2 to z=1.92. However, the light ray generated by yt gives spans from z=2 to z=1.97.
As I look into this, I found that the function comoving_radial_distance in cosmology.py might return the wrong unit. I think the return value should be in comoving units, instead of physical units. To see this directly, I made the following change in cosmology_splice.py:
(root) ~/yt-conda/src/yt-hg/yt/analysis_modules/cosmological_observation $diff cosmology_splice.py cosmology_splice.py0 373d372 < print target_distance, distance2, z2
And it shows:
80.0 Mpccm/h 4.65320708035e+26 cm 1.97387969592 dimensionless 80.0 Mpccm/h 1.19455177073e+26 cm 1.97343368065 dimensionless 80.0 Mpccm/h 1.21507519065e+26 cm 1.97342592996 dimensionless
Then I made the following change in cosmology.py:
(root) ~/yt-conda/src/yt-hg/yt/utilities $diff cosmology.py cosmology.py0 111,112c111,112 < return self.quan((self.hubble_distance() * < trapzint(self.inverse_expansion_factor, z_i,
return (self.hubble_distance() *
trapzint(self.inverse_expansion_factor, z_i,
z_f)).in_base(self.unit_system)
Then I get: 80.0 Mpccm/h 4.65320708035e+26 cmcm 1.92163908776 dimensionless 80.0 Mpccm/h 3.62771515661e+26 cmcm 1.92124702027 dimensionless
With the change of unit from cm to cmcm, the light rays have the right span.
Even though this solves my problem, I am not sure if similar problem still exists. For example, instead of making change in comoving_radial_distance, we might need to change hubble_distance into comoving units. Hopefully someone familiar with yt unit system could check this.
Thanks, Pengfei
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
Hi Britton,
Thanks for your test! In the light ray problem, I think the final redshift should be less than the redshift of the data dump.
My test is like this:
In [5]: from yt.utilities.cosmology import Cosmology ...: cos=Cosmology(hubble_constant = 0.5, ...: omega_matter = 1, ...: omega_lambda = 0, ...: omega_curvature = 0.0) ...: cos.comoving_radial_distance(9.97888260334447,9.9910277956689).to('Mpccm/h')
Out[5]: 0.999963903657 Mpccm/h
Again. as you pointed out, to("Mpc/h") and to("Mpccm/h") would give the same numerical result.
Thanks, Pengfei
On Tue, Aug 9, 2016 at 4:40 AM, Britton Smith brittonsmith@gmail.com wrote:
Hi Pengfei,
Firstly, Nathan, JohnZ, Matt, Cameron, if you could have a look at the last paragraph of this response, I would appreciate it.
Ah, I see now. I actually get a different answer for the redshift interval when I do it by hand. I get 9.994 instead of 9.979. Still, our answers, differ from the value calculated by the code for the reason you describe, essentially that the comoving to proper conversion is different between a cosmology object created from a dataset (where a specific redshift is defined) and one created by hand. In the one created from a dataset, the comoving and proper frame differ by the factor of (1+z), whereas in the one created by hand, they are the same for the reasons I described in my previous email.
This is definitely a problem. We cannot be doing cosmology distance calculations in two different unit systems, one where comoving and proper are the same and one where they are different. As strange as this may sound, the one that is least correct is using the system where the comoving and proper frames are different. That's because something like a comoving radial distance can only be expressed in the comoving frame and has no direct translation to the proper frame via some factor of (1 + z).
The solution to this, as I see it, is the following: The cosmology object that gets associated with a dataset must have its unit_registry overridden to set the comoving and proper frames equal to each other. A better solution would be to even remove the proper frame entirely from the cosmology object's unit_registry, though I don't even know if that's possible. I would appreciate the input of someone with more intimate knowledge of the unit system, like Nathan, Matt, or John Z. As a quick test, I experimented with the case where we create a cosmology object from a dataset by passing a copy of the unit_registry and then modifying the comoving frame. I ran into some difficulties because it looks like, for example, Mpc is define independently of pc, and so I would have to modify both pccm and Mpcm, and perhaps a number of others that I am unaware of. If anyone has ideas on how best to do this, please say so.
In conclusion, I do think there is an issue that needs to be fixed. Input here would be appreciated. Thanks, Pengfei, for identifying this issue!
Britton
On Mon, Aug 8, 2016 at 10:25 PM, Pengfei Chen madcpf@gmail.com wrote:
Hi Britton,
Thanks very much for your test! The error happens when I make light rays from one single dataset. I tried to reproduce the error with the dataset under src/yt-hg/tests. The code I used is here:
http://paste.yt-project.org/show/6751/ After running the code, I did
h5ls -d lightray.h5/grid/redshift and found the last redshift was 9.98999107782409. However, the correct value should be ~ 9.979.
When I call LightRay this way, it will call the _deltaz_forward function in cosmology_splice.py, which will call the comoving_radial_distance in cosmology.py:
distance2 = self.cosmology.comoving_radial_distance(z2, z)
, where distance2 is not in comoving units. Then it will calculate z2:
z2 = ((target_distance - distance2) / m.in_units("Mpccm /
h")) + z2
. I think the subtraction here is the problem, since it tries to convert distance2 to comoving units, even though the numerical value of distance2 is already comoving. So a wrong factor of (1+z) is multiplied.
Any comment is appreciated.
Thanks again for your help! Pengfei
On Mon, Aug 8, 2016 at 7:10 AM, Britton Smith brittonsmith@gmail.com wrote:
Hi Pengfei,
I'm not able to reproduce what you're seeing, but before I go into that, I can explain a few things about how the cosmology calculator works. In yt, we are only able to define a comoving unit system in relation to a proper unit system at a specific redshift. We are unable to define a comoving unit system that is the base unit system. This is a problem for the cosmology calculator, which only calculates comoving quantities. Therefore, in this case, we are forced to set the comoving and proper unit systems to be the same. You can see this for yourself by running the following script:
from yt.utilities.cosmology import Cosmology co = Cosmology() print (co.comoving_radial_distance(1.92, 2).to("Mpc/h")) print (co.comoving_radial_distance(1.92, 2).to("Mpccm/h"))
Both will be the same answer.
Now, back to the issue you're seeing. I ran the following script: http://paste.yt-project.org/show/6741/ to calculate the redshifts dumps needed for a 80 Mpc/h box to go from z = 2 to z = 1.92 and got the following result: CosmologyOutputRedshift[0] = 2.000 CosmologyOutputRedshift[1] = 1.926
This is pretty close to what you expected. Additionally, if I run the following script, http://paste.yt-project.org/show/6743/ with the above redshifts added to the parameter file, I get that dz_max (the change in redshift when moving 80 Mpccm/h at z = 0) of 0.0746521607364, which also seems to agree with what you are finding.
Can you share with me the simulation parameter file and script that you're using to make your light rays? I can take a look at that and see if I still get what I expect.
Britton
On Mon, Aug 8, 2016 at 6:13 AM, Pengfei Chen madcpf@gmail.com wrote:
Hi all,
I was trying to generate light rays using yt 3.4-dev, but found that the light rays have wrong redshift intervals. For example, my simulation box is 80Mpccm/h per side, and I expect a light ray generated in a data dump at z=2 with length =1 simulation unit (i.e. 80Mpccm/h) would have a redshift interval from z=2 to z=1.92. However, the light ray generated by yt gives spans from z=2 to z=1.97.
As I look into this, I found that the function comoving_radial_distance in cosmology.py might return the wrong unit. I think the return value should be in comoving units, instead of physical units. To see this directly, I made the following change in cosmology_splice.py:
(root) ~/yt-conda/src/yt-hg/yt/analysis_modules/cosmological_observation $diff cosmology_splice.py cosmology_splice.py0 373d372 < print target_distance, distance2, z2
And it shows:
80.0 Mpccm/h 4.65320708035e+26 cm 1.97387969592 dimensionless 80.0 Mpccm/h 1.19455177073e+26 cm 1.97343368065 dimensionless 80.0 Mpccm/h 1.21507519065e+26 cm 1.97342592996 dimensionless
Then I made the following change in cosmology.py:
(root) ~/yt-conda/src/yt-hg/yt/utilities $diff cosmology.py cosmology.py0 111,112c111,112 < return self.quan((self.hubble_distance() * < trapzint(self.inverse_expansion_factor, z_i,
return (self.hubble_distance() *
trapzint(self.inverse_expansion_factor, z_i,
z_f)).in_base(self.unit_system)
Then I get: 80.0 Mpccm/h 4.65320708035e+26 cmcm 1.92163908776 dimensionless 80.0 Mpccm/h 3.62771515661e+26 cmcm 1.92124702027 dimensionless
With the change of unit from cm to cmcm, the light rays have the right span.
Even though this solves my problem, I am not sure if similar problem still exists. For example, instead of making change in comoving_radial_distance, we might need to change hubble_distance into comoving units. Hopefully someone familiar with yt unit system could check this.
Thanks, Pengfei
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org
yt-users mailing list yt-users@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-users-spacepope.org