bug in SFRD calculation
There seems to be a bug in the star formation rate density calculation routine. The volume is only the real 'proper volume' at the current redshift and not at any other redshifts, because it uses the current redshift in the calculation of proper volume in all redshifts. #--------------------------------- from yt.mods import * from yt.analysis_modules.star_analysis.api import * import pylab pf=load('HD1200/HD1200') dd=pf.h.all_data() sfr = StarFormationRate(pf, data_source=dd, bins=10) sfr.write_out(name="StarFormationRate.out") time, lookback, z, Msolyr, MsolyrMpc3, Msol, cumMsol=na.loadtxt('StarFormationRate.out',usecols=(0,1,2,3,4,5,6),skiprows=2,unpack=True) cMpc3,=pylab.plot(z,Msolyr/(pf['CosmologyComovingBoxSize']/pf['CosmologyHubbleConstantNow'])**3) currentz,=pylab.plot(z,MsolyrMpc3/(1+pf.current_redshift)**3+1e-5) Expect,=pylab.plot(z,MsolyrMpc3/(1+z)**3) pylab.legend([cMpc3,currentz,Expect], ['[s$^{-1}$cMpc$^{-3}$]', 'use only current z to fit', 'Expected to fit']) pylab.savefig('SFRD.eps') #--------------------------------- in the file: yt/analysis_modules/star_analysis/sfr_spectrum.py line 151: self.Msol_yr_vol.append(self.mass_bins[i] / \ (self.time_bins_dt[i] * tc / YEAR) / vol) The vol is a fixed number instead of varying by the redshift, even though the numerator in solar masses per year is from different redshifts. So I think a couple lines above it should create an empty array of vol vol = [] In the loop for i, time... calculate the right volume vol[i] Then instead of dividing by the constant vol, divide by vol[i] corresponding to the right redshift. self.Msol_yr_vol.append(self.mass_bins[i] / \ (self.time_bins_dt[i] * tc / YEAR) / vol[i]) From G.S.
Hi Geoffrey, I'm not very familiar with the star formation rate calculator, but if you think there is a bug I'd encourage you to submit a pull request with a fix. -Nathan On Thu, Apr 3, 2014 at 6:37 PM, Geoffrey So <gsiisg@gmail.com> wrote:
There seems to be a bug in the star formation rate density calculation routine. The volume is only the real 'proper volume' at the current redshift and not at any other redshifts, because it uses the current redshift in the calculation of proper volume in all redshifts.
#--------------------------------- from yt.mods import * from yt.analysis_modules.star_analysis.api import * import pylab
pf=load('HD1200/HD1200') dd=pf.h.all_data() sfr = StarFormationRate(pf, data_source=dd, bins=10) sfr.write_out(name="StarFormationRate.out") time, lookback, z, Msolyr, MsolyrMpc3, Msol, cumMsol=na.loadtxt('StarFormationRate.out',usecols=(0,1,2,3,4,5,6),skiprows=2,unpack=True)
cMpc3,=pylab.plot(z,Msolyr/(pf['CosmologyComovingBoxSize']/pf['CosmologyHubbleConstantNow'])**3) currentz,=pylab.plot(z,MsolyrMpc3/(1+pf.current_redshift)**3+1e-5) Expect,=pylab.plot(z,MsolyrMpc3/(1+z)**3)
pylab.legend([cMpc3,currentz,Expect], ['[s$^{-1}$cMpc$^{-3}$]', 'use only current z to fit', 'Expected to fit'])
pylab.savefig('SFRD.eps') #---------------------------------
in the file: yt/analysis_modules/star_analysis/sfr_spectrum.py line 151: self.Msol_yr_vol.append(self.mass_bins[i] / \ (self.time_bins_dt[i] * tc / YEAR) / vol)
The vol is a fixed number instead of varying by the redshift, even though the numerator in solar masses per year is from different redshifts.
So I think a couple lines above it should create an empty array of vol vol = [] In the loop for i, time... calculate the right volume vol[i] Then instead of dividing by the constant vol, divide by vol[i] corresponding to the right redshift. self.Msol_yr_vol.append(self.mass_bins[i] / \ (self.time_bins_dt[i] * tc / YEAR) / vol[i])
From G.S.
_______________________________________________ yt-dev mailing list yt-dev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org
Yes, although the star formation rate in cosmology is usually calculated in a comoving volume, not proper (see, e.g., http://arxiv.org/pdf/1403.0007v2.pdf). Greg On Apr 3, 2014, at 11:15 PM, Nathan Goldbaum <nathan12343@gmail.com> wrote:
Hi Geoffrey,
I'm not very familiar with the star formation rate calculator, but if you think there is a bug I'd encourage you to submit a pull request with a fix.
-Nathan
On Thu, Apr 3, 2014 at 6:37 PM, Geoffrey So <gsiisg@gmail.com> wrote: There seems to be a bug in the star formation rate density calculation routine. The volume is only the real 'proper volume' at the current redshift and not at any other redshifts, because it uses the current redshift in the calculation of proper volume in all redshifts.
#--------------------------------- from yt.mods import * from yt.analysis_modules.star_analysis.api import * import pylab
pf=load('HD1200/HD1200') dd=pf.h.all_data() sfr = StarFormationRate(pf, data_source=dd, bins=10) sfr.write_out(name="StarFormationRate.out") time, lookback, z, Msolyr, MsolyrMpc3, Msol, cumMsol=na.loadtxt('StarFormationRate.out',usecols=(0,1,2,3,4,5,6),skiprows=2,unpack=True)
cMpc3,=pylab.plot(z,Msolyr/(pf['CosmologyComovingBoxSize']/pf['CosmologyHubbleConstantNow'])**3) currentz,=pylab.plot(z,MsolyrMpc3/(1+pf.current_redshift)**3+1e-5) Expect,=pylab.plot(z,MsolyrMpc3/(1+z)**3)
pylab.legend([cMpc3,currentz,Expect], ['[s$^{-1}$cMpc$^{-3}$]', 'use only current z to fit', 'Expected to fit'])
pylab.savefig('SFRD.eps') #---------------------------------
in the file: yt/analysis_modules/star_analysis/sfr_spectrum.py line 151: self.Msol_yr_vol.append(self.mass_bins[i] / \ (self.time_bins_dt[i] * tc / YEAR) / vol)
The vol is a fixed number instead of varying by the redshift, even though the numerator in solar masses per year is from different redshifts.
So I think a couple lines above it should create an empty array of vol vol = [] In the loop for i, time… calculate the right volume vol[i] Then instead of dividing by the constant vol, divide by vol[i] corresponding to the right redshift. self.Msol_yr_vol.append(self.mass_bins[i] / \ (self.time_bins_dt[i] * tc / YEAR) / vol[i])
From G.S.
_______________________________________________ yt-dev mailing list yt-dev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org
_______________________________________________ yt-dev mailing list yt-dev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org
---------------------------------- Greg Bryan Associate Professor of Astronomy Department of Astronomy, Columbia University 1325 Pupin Physics Laboratories, Mail Code 5246 550 West 120th Street New York, New York 10027 (USA) email: gbryan@astro.columbia.edu phone: +1(212)854-6837 fax : +1(212)854-8121
Hi Greg, My original script has them converted from proper to comoving and that's when I noticed the bug :) Nathan, I was trying to convey the idea, wasn't sure if what I described is the proper fix, I'll issue the PR if it tests out correctly. From G.S. Sent from my iPhone
On Apr 4, 2014, at 6:15 AM, Greg Bryan <gbryan@astro.columbia.edu> wrote:
Yes, although the star formation rate in cosmology is usually calculated in a comoving volume, not proper (see, e.g., http://arxiv.org/pdf/1403.0007v2.pdf).
Greg
On Apr 3, 2014, at 11:15 PM, Nathan Goldbaum <nathan12343@gmail.com> wrote:
Hi Geoffrey,
I'm not very familiar with the star formation rate calculator, but if you think there is a bug I'd encourage you to submit a pull request with a fix.
-Nathan
On Thu, Apr 3, 2014 at 6:37 PM, Geoffrey So <gsiisg@gmail.com> wrote: There seems to be a bug in the star formation rate density calculation routine. The volume is only the real 'proper volume' at the current redshift and not at any other redshifts, because it uses the current redshift in the calculation of proper volume in all redshifts.
#--------------------------------- from yt.mods import * from yt.analysis_modules.star_analysis.api import * import pylab
pf=load('HD1200/HD1200') dd=pf.h.all_data() sfr = StarFormationRate(pf, data_source=dd, bins=10) sfr.write_out(name="StarFormationRate.out") time, lookback, z, Msolyr, MsolyrMpc3, Msol, cumMsol=na.loadtxt('StarFormationRate.out',usecols=(0,1,2,3,4,5,6),skiprows=2,unpack=True)
cMpc3,=pylab.plot(z,Msolyr/(pf['CosmologyComovingBoxSize']/pf['CosmologyHubbleConstantNow'])**3) currentz,=pylab.plot(z,MsolyrMpc3/(1+pf.current_redshift)**3+1e-5) Expect,=pylab.plot(z,MsolyrMpc3/(1+z)**3)
pylab.legend([cMpc3,currentz,Expect], ['[s$^{-1}$cMpc$^{-3}$]', 'use only current z to fit', 'Expected to fit'])
pylab.savefig('SFRD.eps') #---------------------------------
in the file: yt/analysis_modules/star_analysis/sfr_spectrum.py line 151: self.Msol_yr_vol.append(self.mass_bins[i] / \ (self.time_bins_dt[i] * tc / YEAR) / vol)
The vol is a fixed number instead of varying by the redshift, even though the numerator in solar masses per year is from different redshifts.
So I think a couple lines above it should create an empty array of vol vol = [] In the loop for i, time… calculate the right volume vol[i] Then instead of dividing by the constant vol, divide by vol[i] corresponding to the right redshift. self.Msol_yr_vol.append(self.mass_bins[i] / \ (self.time_bins_dt[i] * tc / YEAR) / vol[i])
From G.S.
_______________________________________________ yt-dev mailing list yt-dev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org
_______________________________________________ yt-dev mailing list yt-dev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org
---------------------------------- Greg Bryan Associate Professor of Astronomy Department of Astronomy, Columbia University 1325 Pupin Physics Laboratories, Mail Code 5246 550 West 120th Street New York, New York 10027 (USA)
email: gbryan@astro.columbia.edu phone: +1(212)854-6837 fax : +1(212)854-8121
_______________________________________________ yt-dev mailing list yt-dev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org
Hi, I've issued a pull request to yt_analysis, the change diverting from the norm for yt to output in proper units, but as Greg pointed out in the case of star formation rate density everyone uses comoving Mpc, so I think the change will be useful. The original quantity doesn't make much sense, it was dividing Msol/yr at all redshifts by the proper volume at the last redshift. The script should show identical curves offset by a tiny bit. From G.S. #------------ from yt.mods import * from yt.analysis_modules.star_analysis.api import * import pylab pf=load('HD1200/HD1200') dd=pf.h.all_data() sfr = StarFormationRate(pf, data_source=dd, bins=10) sfr.write_out(name="StarFormationRate.out") time, lookback, z, Msolyr, MsolyrMpc3, Msol, cumMsol=na.loadtxt('StarFormationRate.out',usecols=(0,1,2,3,4,5,6),skiprows=2,unpack=True) cMpc3,=pylab.plot(z,Msolyr/(pf['CosmologyComovingBoxSize']/pf['CosmologyHubbleConstantNow'])**3) New,=pylab.plot(z,MsolyrMpc3+3e-5) pylab.legend([cMpc3,NewMpc3], ['[s$^{-1}$cMpc$^{-3}$]', 'New']) pylab.savefig('SFRD.eps') #------------ On Fri, Apr 4, 2014 at 9:36 AM, Geoffrey So <gsiisg@gmail.com> wrote:
Hi Greg,
My original script has them converted from proper to comoving and that's when I noticed the bug :)
Nathan,
I was trying to convey the idea, wasn't sure if what I described is the proper fix, I'll issue the PR if it tests out correctly.
From G.S.
Sent from my iPhone
On Apr 4, 2014, at 6:15 AM, Greg Bryan <gbryan@astro.columbia.edu> wrote:
Yes, although the star formation rate in cosmology is usually calculated in a comoving volume, not proper (see, e.g., http://arxiv.org/pdf/1403.0007v2.pdf).
Greg
On Apr 3, 2014, at 11:15 PM, Nathan Goldbaum <nathan12343@gmail.com> wrote:
Hi Geoffrey,
I'm not very familiar with the star formation rate calculator, but if you think there is a bug I'd encourage you to submit a pull request with a fix.
-Nathan
On Thu, Apr 3, 2014 at 6:37 PM, Geoffrey So <gsiisg@gmail.com> wrote:
There seems to be a bug in the star formation rate density calculation routine. The volume is only the real 'proper volume' at the current redshift and not at any other redshifts, because it uses the current redshift in the calculation of proper volume in all redshifts.
#--------------------------------- from yt.mods import * from yt.analysis_modules.star_analysis.api import * import pylab
pf=load('HD1200/HD1200') dd=pf.h.all_data() sfr = StarFormationRate(pf, data_source=dd, bins=10) sfr.write_out(name="StarFormationRate.out") time, lookback, z, Msolyr, MsolyrMpc3, Msol, cumMsol=na.loadtxt('StarFormationRate.out',usecols=(0,1,2,3,4,5,6),skiprows=2,unpack=True)
cMpc3,=pylab.plot(z,Msolyr/(pf['CosmologyComovingBoxSize']/pf['CosmologyHubbleConstantNow'])**3) currentz,=pylab.plot(z,MsolyrMpc3/(1+pf.current_redshift)**3+1e-5) Expect,=pylab.plot(z,MsolyrMpc3/(1+z)**3)
pylab.legend([cMpc3,currentz,Expect], ['[s$^{-1}$cMpc$^{-3}$]', 'use only current z to fit', 'Expected to fit'])
pylab.savefig('SFRD.eps') #---------------------------------
in the file: yt/analysis_modules/star_analysis/sfr_spectrum.py line 151: self.Msol_yr_vol.append(self.mass_bins[i] / \ (self.time_bins_dt[i] * tc / YEAR) / vol)
The vol is a fixed number instead of varying by the redshift, even though the numerator in solar masses per year is from different redshifts.
So I think a couple lines above it should create an empty array of vol vol = [] In the loop for i, time... calculate the right volume vol[i] Then instead of dividing by the constant vol, divide by vol[i] corresponding to the right redshift. self.Msol_yr_vol.append(self.mass_bins[i] / \ (self.time_bins_dt[i] * tc / YEAR) / vol[i])
From G.S.
_______________________________________________ yt-dev mailing list yt-dev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org
_______________________________________________ yt-dev mailing list yt-dev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org
---------------------------------- Greg Bryan Associate Professor of Astronomy Department of Astronomy, Columbia University 1325 Pupin Physics Laboratories, Mail Code 5246 550 West 120th Street New York, New York 10027 (USA)
email: gbryan@astro.columbia.edu phone: +1(212)854-6837 fax : +1(212)854-8121
_______________________________________________ yt-dev mailing list yt-dev@lists.spacepope.org http://lists.spacepope.org/listinfo.cgi/yt-dev-spacepope.org
participants (3)
-
Geoffrey So -
Greg Bryan -
Nathan Goldbaum