astropy / halotools

Python package for studying large scale structure, cosmology, and galaxy evolution using N-body simulations and halo models
http://halotools.rtfd.org
99 stars 63 forks source link

case use: DeltaSigma when galaxy-matter cross tpcf is negative #216

Closed duncandc closed 8 years ago

duncandc commented 8 years ago

It is possible that the galaxy matter cross correlation, xi_gm, takes negative values.

Currently we integrate the log(xi_gm) when calculating the galaxy-galaxy lensing signal DeltaSigma.

Right now, we throw an error if xi_gm<0. Is there something better to do?

aphearin commented 8 years ago
  1. Integrate in linear space
  2. Perform the integral in log-space but do it piece-wise
  3. Shift the values to be positive before logarithmic integration, then shift them back.

@surhudm @vdbosch69 - other suggestions?

surhudm commented 8 years ago

Can you point me to the code where this is done? I do not understand why you are integrating log(xi_gm) at all.

On Thu, Nov 19, 2015 at 11:28 PM, Andrew Hearin notifications@github.com wrote:

  1. Integrate in linear space
  2. Perform the integral in log-space but do it piece-wise
  3. Shift the values to be positive before logarithmic integration, then shift them back.

@surhudm https://github.com/surhudm @vdbosch69 https://github.com/vdbosch69 - other suggestions?

— Reply to this email directly or view it on GitHub https://github.com/astropy/halotools/issues/216#issuecomment-158072832.

vdbosch69 commented 8 years ago

haven’t thought about this much, but perhaps you can integrate 1+xi_gm(r) using logarithmic integration over r, and then subtract the 1+ component (this should work unless you integrate out to infinity).

F

On Nov 19, 2015, at 9:28 AM, Andrew Hearin notifications@github.com wrote:

Integrate in linear space Perform the integral in log-space but do it piece-wise Shift the values to be positive before logarithmic integration, then shift them back. @surhudm https://urldefense.proofpoint.com/v2/url?u=https-3A__github.com_surhudm&d=AwMCaQ&c=-dg2m7zWuuDZ0MUcV7Sdqw&r=akBy6VKVSXuVF7qPhX_uKgCTzeIo3SdhrMQi2H7IAuI&m=cOGj06ejoykQ8Vo8W0Qy3kAg1d5USr89LGjhq2o7mhs&s=PVwOhRHjZiHEF-iy0IntEf1PEaFrzXZALeu-cg5Uc6Y&e= @vdbosch69 https://urldefense.proofpoint.com/v2/url?u=https-3A__github.com_vdbosch69&d=AwMCaQ&c=-dg2m7zWuuDZ0MUcV7Sdqw&r=akBy6VKVSXuVF7qPhX_uKgCTzeIo3SdhrMQi2H7IAuI&m=cOGj06ejoykQ8Vo8W0Qy3kAg1d5USr89LGjhq2o7mhs&s=ZnfNxZRVcRlKAfXeprs--9wIKIOPegb0kPh5AHY9ySs&e= - other suggestions?

— Reply to this email directly or view it on GitHub https://urldefense.proofpoint.com/v2/url?u=https-3A__github.com_astropy_halotools_issues_216-23issuecomment-2D158072832&d=AwMCaQ&c=-dg2m7zWuuDZ0MUcV7Sdqw&r=akBy6VKVSXuVF7qPhX_uKgCTzeIo3SdhrMQi2H7IAuI&m=cOGj06ejoykQ8Vo8W0Qy3kAg1d5USr89LGjhq2o7mhs&s=kAHITOR6NKUic9dm0g4BmG8FfUBkRAiLtfamAG72Mc0&e=.

vdbosch69 commented 8 years ago

I think he means he is integrating xi_gm but using logarithmic integration intervals in r. Is that correct AH?

On Nov 19, 2015, at 9:38 AM, Surhud More notifications@github.com wrote:

Can you point me to the code where this is done? I do not understand why you are integrating log(xi_gm) at all.

On Thu, Nov 19, 2015 at 11:28 PM, Andrew Hearin notifications@github.com wrote:

  1. Integrate in linear space
  2. Perform the integral in log-space but do it piece-wise
  3. Shift the values to be positive before logarithmic integration, then shift them back.

@surhudm https://github.com/surhudm @vdbosch69 https://github.com/vdbosch69 - other suggestions?

— Reply to this email directly or view it on GitHub https://github.com/astropy/halotools/issues/216#issuecomment-158072832.

— Reply to this email directly or view it on GitHub https://urldefense.proofpoint.com/v2/url?u=https-3A__github.com_astropy_halotools_issues_216-23issuecomment-2D158075208&d=AwMFaQ&c=-dg2m7zWuuDZ0MUcV7Sdqw&r=akBy6VKVSXuVF7qPhX_uKgCTzeIo3SdhrMQi2H7IAuI&m=RHbhwF5Iiptoqe4pLlc-l5tkmyWbt8HcQ4zfpm6c538&s=X0DUsiGtn7FTqsiMH_DhqdUdAmRKqHhZm1R0R2Pq5gE&e=.

aphearin commented 8 years ago

This code is not in the master branch - @duncandc has been developing it in a branch on his fork so we'll have to wait for him to point us to the specific branch. However, now that you point it out I see that I didn't understand this Issue before commenting. I also don't know where taking the log of xi_gm comes into play when computing DS. Let's wait for @duncandc to clarify.

surhudm commented 8 years ago

Ok I found the branch then, @duncandc is basically interpolating log(xigm) vs log(r), and that will pose a problem if you have negative xigm. The solution is simple, interpolate log(1+xigm) vs log(r), this is exactly what we do in the codes that @vdbosch69 and I developed and what I mentioned in my comment on #218.

On Thu, Nov 19, 2015 at 11:57 PM, Andrew Hearin notifications@github.com wrote:

This code is not in the master branch - @duncandc https://github.com/duncandc has been developing it in a branch on his fork so we'll have to wait for him to point us to the specific branch. However, now that you point it out I see that I didn't understand this Issue before commenting. I also don't know where taking the log of xi_gm comes into play when computing DS. Let's wait for @duncandc https://github.com/duncandc to clarify.

— Reply to this email directly or view it on GitHub https://github.com/astropy/halotools/issues/216#issuecomment-158081230.

duncandc commented 8 years ago

https://github.com/duncandc/halotools/blob/mock_obs_devel/halotools/mock_observables/delta_sigma.py

aphearin commented 8 years ago

Right, so your integrand is a spline of xi(r) and you build the spline in log-log space. Are you sure log-log interpolation is necessary? Do you run into numerical problems if you build the table in log-linear space?

aphearin commented 8 years ago

Brief aside: this code is exceptionally clearly written and easy to understand. Really nice work.

duncandc commented 8 years ago

I am going to fit log(1+xi). Any other issues, we can open a new issue thread.

manodeep commented 8 years ago

Seems like this issue is different from the one I have seen. Does the tpcf use the binning in dec/ra ?

duncandc commented 8 years ago

@manodeep what issue are you referring to? This uses the real space tpcf function. binning is done in r. the pair counter uses cells.

manodeep commented 8 years ago

So the inputs positions are (xyz) - rather than ra/dec/cz? For some of Jen's weird correlation functions on mocks, I got negative xi which went away once I removed the binning in my codes. In all of those cases, pimax was effectively infinity..

On a different note, you can improve performance by switching the interpolation from r to r^2, i.e.,

set up the initial interpolation as:

 xi = InterpolatedUnivariateSpline(rbin_centers**2, np.log10(xi+1.0), ext=0)

 def f(pi,rp):
     r2 = rp**2+pi**2
     return mean_rho*(1.0+(10.0**xi(r2)-1.0))

and then you can avoid the sqrt inside f().