First import all the modules such as healpy and astropy needed for analyzing the structure
In [2]:
import healpix_util as hu
import astropy as ap
import numpy as np
from astropy.io import fits
from astropy.table import Table
import astropy.io.ascii as ascii
from astropy.io import fits
from astropy.constants import c
import matplotlib.pyplot as plt
import math as m
from math import pi
#from scipy.constants import c
import scipy.special as sp
from astroML.decorators import pickle_results
from scipy import integrate
import warnings
from sklearn.neighbors import BallTree
import pickle
import multiprocessing as mp
import time
from aptestmetricdt import *
from aptestmetricdz import *
from scipy.spatial import distance as d
from apcat import *
from progressbar import *
from tqdm import *
from functools import partial
import pymangle
from scipy.optimize import curve_fit
#from astroML.datasets import fetch_sdss_specgals
#from astroML.correlation import bootstrap_two_point_angular
%matplotlib inline
Read the data file (taken from http://cosmo.nyu.edu/~eak306/SDSS-LRG.html ) converted to ascii with comoving distance etc. in V01 reading from pkl files for faster read
In [ ]:
In [3]:
with open('dd2ddr72v06cdist200kopt2.pkl') as f:
dd2d=pickle.load(f)
In [4]:
with open('rr2ddr72v06cdist200kopt2.pkl') as f:
rr2d=pickle.load(f)
In [5]:
with open('dr2ddr72v06cdist200kopt2.pkl') as f:
dr2d=pickle.load(f)
In [13]:
dd2d
Out[13]:
In [14]:
rr2d
Out[14]:
In [15]:
dr2d
Out[15]:
In [16]:
corrells=(4.0*dd2d+rr2d*1.0-4.0*dr2d)/rr2d
In [17]:
corrells
Out[17]:
In [40]:
sigp=(1.0+corrells)/dd2d
In [52]:
sigp
Out[52]:
In [55]:
ans3=curve_fit(xi0,dzbin[1:],corr0,bounds=(0, [3., 2., 1.]),maxfev=10000,sigma=sigp[0][1:])
In [56]:
ans3
Out[56]:
In [58]:
err3=np.sqrt(np.diag(ans3[1]))
err3
Out[58]:
In [59]:
plt.plot(dzbin[1:],corr0,'bo-')
plt.plot(dzbin[1:],xi0(dzbin[1:],1.72940481, 0.75864624, 0.18241307),'ro-')
plt.plot(dzbin[1:],xi0(dzbin[1:],1.88029658, 0.57632231, 0.57244037),'go-')
Out[59]:
In [ ]:
In [54]:
sigp.shape
Out[54]:
In [44]:
import scipy
In [47]:
scipy.__version__
Out[47]:
In [8]:
plt.contour(corrells)
Out[8]:
In [9]:
dzbin=zdthbin=np.arange(0.001,0.0201,0.001)
In [10]:
dzbin
Out[10]:
In [11]:
plt.contour(dzbin,zdthbin,corrells)
Out[11]:
In [12]:
ztotal = np.concatenate([np.fliplr(corrells),corrells],axis=1)
ztotal = np.concatenate([np.flipud(ztotal),ztotal],axis=0)
xtotal = np.concatenate([-dzbin[::-1],zdthbin],axis=0)
ytotal = np.concatenate([-zdthbin[::-1],zdthbin],axis=0)
plt.figure()
plt.contour(xtotal,ytotal,ztotal)
#plt.savefig("../plots/anisotropic2pcfdzzdth.pdf")
plt.show()
In [18]:
cor0=np.array([[ 3.21829653e+00, 5.15924709e+00, 2.79014482e+00,
1.61622821e+00, 1.08643298e+00, 7.84456307e-01,
5.49128339e-01, 4.85461602e-01, 3.47976294e-01,
2.64233207e-01, 2.64375860e-01, 1.98222947e-01,
1.24068892e-01, 1.18280361e-01, 7.35890402e-02,
7.78017101e-02, 8.63556000e-02, 6.77324604e-02,
5.80771512e-02, 6.67772329e-02],
[ 8.10972346e+00, 3.41727694e+00, 2.04130609e+00,
1.34601758e+00, 9.73548486e-01, 7.11285425e-01,
5.49041632e-01, 4.32106559e-01, 3.39539082e-01,
2.96327191e-01, 2.24821404e-01, 2.02357985e-01,
1.29825807e-01, 1.19599598e-01, 9.34292995e-02,
1.24206445e-01, 1.08218439e-01, 9.93404422e-02,
7.22001764e-02, 6.65998734e-02],
[ 4.84014210e+00, 1.95147392e+00, 1.33814659e+00,
1.06878030e+00, 8.15801446e-01, 6.21011541e-01,
5.21177836e-01, 3.97793227e-01, 2.95925863e-01,
2.48257618e-01, 2.20741774e-01, 1.79282525e-01,
1.47969417e-01, 1.16869617e-01, 1.38469511e-01,
9.14332248e-02, 1.10308692e-01, 8.89953816e-02,
9.02822578e-02, 6.54145982e-02],
[ 2.91616209e+00, 1.18162932e+00, 8.81038212e-01,
7.25515886e-01, 6.12496184e-01, 4.83955635e-01,
3.90963392e-01, 3.10118435e-01, 2.43897553e-01,
2.46864076e-01, 1.63464464e-01, 1.49259531e-01,
1.38636227e-01, 1.43624410e-01, 1.14456475e-01,
9.38243441e-02, 9.56588789e-02, 8.28049616e-02,
7.09190853e-02, 6.38867889e-02],
[ 1.85621521e+00, 7.88299156e-01, 6.62597522e-01,
4.70473369e-01, 4.65349607e-01, 4.26253875e-01,
3.13391151e-01, 2.41903060e-01, 2.37370559e-01,
1.94981949e-01, 1.96513545e-01, 1.41195546e-01,
1.11824162e-01, 1.12188597e-01, 8.54642987e-02,
1.05656979e-01, 9.28853449e-02, 7.73981062e-02,
8.68920377e-02, 5.88907947e-02],
[ 1.10414828e+00, 5.30405405e-01, 4.26055046e-01,
4.02095828e-01, 3.35974148e-01, 3.08711465e-01,
2.78604600e-01, 2.66722408e-01, 2.17977271e-01,
1.79992811e-01, 1.58517890e-01, 1.39719525e-01,
1.33705743e-01, 1.15663504e-01, 1.05410317e-01,
9.18325416e-02, 6.66625335e-02, 7.44469173e-02,
7.50256072e-02, 8.50809980e-02],
[ 8.14265199e-01, 3.80754381e-01, 3.17623510e-01,
3.42954605e-01, 2.72279045e-01, 2.06386293e-01,
2.13127037e-01, 1.89496097e-01, 1.57882637e-01,
1.70882808e-01, 1.33449536e-01, 1.17673883e-01,
1.18497168e-01, 1.14211115e-01, 9.38498823e-02,
6.84812702e-02, 9.53366424e-02, 7.04954655e-02,
8.16902195e-02, 6.56066553e-02],
[ 4.80743090e-01, 2.45037037e-01, 2.60309278e-01,
2.04044974e-01, 1.75123762e-01, 2.01413732e-01,
1.95487058e-01, 1.49848505e-01, 1.50522729e-01,
1.36128781e-01, 9.08532262e-02, 1.14173825e-01,
1.12524046e-01, 8.51803316e-02, 9.41305192e-02,
7.09896024e-02, 7.74198580e-02, 8.52549311e-02,
6.40647391e-02, 9.58478019e-02],
[ 4.54304636e-01, 1.77953234e-01, 1.78434847e-01,
1.87964975e-01, 1.85141885e-01, 1.52230290e-01,
1.70238184e-01, 1.15514362e-01, 1.20041543e-01,
1.04337324e-01, 1.25766142e-01, 7.86827056e-02,
7.36177824e-02, 7.03130923e-02, 8.58850152e-02,
6.93294910e-02, 7.42290505e-02, 6.69872178e-02,
6.81848627e-02, 5.86442644e-02],
[ 1.95520299e-01, 2.35311573e-01, 1.07520891e-01,
1.11604440e-01, 1.40897821e-01, 1.28361062e-01,
1.01710104e-01, 1.29342931e-01, 1.08348645e-01,
7.25599875e-02, 1.08945717e-01, 7.11690813e-02,
8.45656099e-02, 6.93269304e-02, 7.87809903e-02,
5.01656571e-02, 8.17190145e-02, 5.72408731e-02,
4.29866850e-02, 6.23349324e-02],
[ 1.30824373e-01, 5.00863558e-02, 1.33339502e-01,
1.11840783e-01, 9.07466765e-02, 1.18902573e-01,
1.04517469e-01, 8.26933852e-02, 8.98120615e-02,
9.75724847e-02, 7.31578187e-02, 7.77130663e-02,
8.66212047e-02, 5.55312137e-02, 9.08994184e-02,
7.74623204e-02, 4.57442289e-02, 3.14908676e-02,
3.57278902e-02, 4.26337493e-02],
[ 2.62302483e-01, 8.08065261e-02, 1.16038433e-01,
8.35241816e-02, 1.06246043e-01, 4.75735741e-02,
1.13633826e-01, 8.14940577e-02, 8.84401522e-02,
7.64548460e-02, 1.02468572e-01, 4.04988572e-02,
8.40946548e-02, 4.37316666e-02, 7.96138639e-02,
4.30466930e-02, 4.18495198e-02, 5.96056508e-02,
4.47325492e-02, 5.11649403e-02],
[ 5.96219098e-02, 6.01728201e-02, 1.16514691e-01,
9.45690354e-02, 9.24812823e-02, 6.84446020e-02,
7.66094015e-02, 7.68422090e-02, 4.21049547e-02,
4.31055013e-02, 5.18608023e-02, 7.79947370e-02,
4.59949809e-02, 3.78655839e-02, 4.13113102e-02,
3.56767745e-02, 4.54018884e-02, 5.05821319e-02,
4.10740017e-02, 4.30719634e-02],
[ 2.30434783e-01, 9.97329983e-02, 9.51249644e-02,
6.11103533e-02, 7.30459801e-02, 5.05657709e-02,
5.51406804e-02, 4.94054267e-02, 8.86476046e-02,
5.96951285e-02, 8.12231708e-02, 5.28811006e-02,
3.44244025e-02, 6.10717780e-02, 5.68130661e-02,
6.83014049e-02, 1.81400469e-02, 5.17377089e-02,
4.82751663e-02, 5.16523700e-02],
[ 0, 9.68351441e-02, 5.28748916e-02,
5.38193257e-02, 2.37574075e-02, 4.50828828e-02,
6.75124604e-02, 5.39701612e-02, 7.58320337e-02,
5.09614875e-02, 3.36042703e-02, 2.21982807e-02,
3.99732676e-02, 3.20129377e-02, 4.11943101e-02,
3.57650362e-02, 5.52296373e-02, 2.13559581e-02,
4.56760576e-02, 3.26732537e-02],
[ 0, 3.72750643e-02, 8.95989381e-02,
3.98218143e-02, 4.49058426e-02, 2.05173952e-02,
4.64911946e-02, 5.95096876e-02, 4.75283474e-02,
4.60800000e-02, 5.98061657e-02, 2.93272254e-02,
4.64577036e-02, 5.76859856e-02, 3.41340995e-02,
4.27923689e-02, 3.43455260e-02, 2.25067612e-02,
3.64251390e-02, 2.85137557e-02],
[ 5.67107750e-02, 5.34184706e-02, 5.34821085e-02,
0, 4.06704991e-02, 0,
9.19978641e-02, 5.18196995e-02, 1.91692335e-02,
5.04086828e-02, 4.25723975e-02, 5.06320587e-02,
6.86611694e-02, 5.20188595e-02, 2.84061116e-02,
2.25367597e-02, 2.75510533e-02, 1.79712460e-02,
2.85947951e-02, 1.42072257e-02],
[ 0, 0, 0,
0, 0, 4.26870066e-02,
5.84179276e-02, 1.87100445e-02, 3.15694528e-02,
3.79542396e-02, 4.59629157e-02, 7.50452080e-02,
3.89305620e-02, 3.28357039e-02, 4.95586242e-02,
2.61048657e-02, 4.08723774e-02, 3.44360462e-02,
3.70191215e-02, 1.71895593e-02],
[ 5.11415525e-02, 0.00000000e+00, 5.24283094e-02,
1.77886292e-02, 1.83382630e-02, 2.58956125e-02,
3.35668083e-02, 3.27660145e-02, 2.92453966e-02,
1.77364408e-02, 7.32768208e-02, 3.79056449e-02,
3.34790782e-02, 5.12010236e-02, 4.72340114e-02,
3.45161849e-02, 3.72760326e-02, 2.75757949e-02,
2.88515160e-02, 1.31145696e-02],
[ 1.33843212e-01, 1.15074980e-01, 2.23581120e-02,
2.61687531e-03, 3.52228949e-03, 1.85335133e-02,
2.60140075e-02, 4.51142895e-02, 3.66334547e-02,
4.30378375e-02, 5.28900457e-02, 5.72485951e-02,
1.63365499e-02, 1.18434604e-02, 2.17499954e-02,
1.06397775e-02, 2.91595759e-02, 2.39908370e-02,
2.47498529e-02, 3.24743642e-02]])
In [19]:
plt.contour(dzbin,zdthbin,cor0)
Out[19]:
In [24]:
corr0=np.array([ 5.15924709e+00, 2.79014482e+00,
1.61622821e+00, 1.08643298e+00, 7.84456307e-01,
5.49128339e-01, 4.85461602e-01, 3.47976294e-01,
2.64233207e-01, 2.64375860e-01, 1.98222947e-01,
1.24068892e-01, 1.18280361e-01, 7.35890402e-02,
7.78017101e-02, 8.63556000e-02, 6.77324604e-02,
5.80771512e-02, 6.67772329e-02])
In [38]:
ans2=curve_fit(xi0,dzbin[1:],corr0,bounds=(0, [3., 2., 1.]),maxfev=10000)
In [39]:
ans2
Out[39]:
In [28]:
ans
Out[28]:
In [30]:
plt.plot(dzbin[1:],corr0,'bo-')
plt.plot(dzbin[1:],xi0(dzbin[1:],1.72940481, 0.75864624, 0.18241307),'ro-')
Out[30]:
In [31]:
perr=np.sqrt(np.diag(ans[1]))
In [32]:
perr
Out[32]:
In [33]:
rd2d=np.array([[ 2205., 6509., 10617., 14888., 19135., 22806., 27350.,
30627., 34655., 38549., 41577., 45369., 49577., 53261.,
57240., 60768., 63806., 67306., 71212., 74766.],
[ 2128., 6383., 10532., 14979., 18853., 23107., 26913.,
30642., 34333., 38065., 41954., 45513., 49499., 53067.,
57135., 60108., 63503., 67014., 70805., 74657.],
[ 2092., 6391., 10654., 14809., 18990., 22834., 26508.,
30568., 34650., 38387., 42084., 45143., 49488., 53330.,
56165., 60097., 63303., 67144., 70479., 74100.],
[ 2247., 6400., 10500., 14779., 18684., 22789., 26594.,
30704., 34519., 38178., 42047., 45556., 48968., 52392.,
56163., 60330., 63307., 66879., 70654., 74341.],
[ 2177., 6230., 10340., 14784., 18631., 22585., 26420.,
30845., 34098., 38171., 41496., 45234., 49094., 52691.,
56408., 59457., 63411., 66856., 70481., 74018.],
[ 2096., 6303., 10624., 14560., 18722., 22623., 26528.,
30239., 34154., 38140., 41459., 45541., 49034., 52193.,
55814., 59624., 63021., 66900., 70365., 73366.],
[ 2116., 6333., 10478., 14634., 18835., 22557., 26436.,
30340., 34052., 37778., 41286., 45360., 48667., 52378.,
56007., 59912., 62522., 66802., 69625., 73669.],
[ 2189., 6326., 10472., 14840., 18712., 22568., 26304.,
30083., 33981., 37972., 41651., 45039., 48446., 52628.,
55612., 59590., 62745., 66145., 70004., 72421.],
[ 2166., 6325., 10496., 14407., 18415., 22493., 26147.,
30209., 33923., 37973., 41162., 45098., 49033., 52131.,
55565., 59018., 62269., 65989., 69622., 72872.],
[ 2100., 6255., 10528., 14459., 18531., 22533., 26487.,
30160., 33683., 37591., 41224., 44711., 48625., 51873.,
55524., 59254., 62227., 65706., 69482., 72907.],
[ 2128., 6457., 10431., 14537., 18634., 22172., 26048.,
29986., 33730., 37514., 40935., 44329., 48276., 52141.,
55179., 58650., 62581., 66163., 69452., 72689.],
[ 2077., 6304., 10282., 14489., 18333., 22617., 25961.,
30042., 33060., 37106., 40536., 44903., 48005., 51465.,
54936., 58906., 62288., 65332., 69156., 72467.],
[ 2096., 6283., 10304., 14390., 18162., 22164., 26022.,
29688., 33532., 37415., 41027., 44497., 48171., 51561.,
55065., 58332., 62011., 65375., 68791., 72176.],
[ 1997., 6202., 10161., 14317., 18073., 22174., 25932.,
29731., 33118., 37008., 40333., 44390., 48149., 50994.,
54619., 57655., 61769., 64973., 68538., 72150.],
[ 2177., 6144., 10177., 14274., 18243., 21932., 25646.,
29484., 32853., 37053., 40587., 44570., 47524., 51278.,
54741., 57904., 60909., 65417., 68491., 72215.],
[ 2144., 6248., 10135., 14215., 18261., 22048., 25676.,
29599., 33066., 36720., 40301., 43999., 47269., 51008.,
54608., 57814., 61531., 65049., 68119., 71577.],
[ 2056., 6089., 10218., 14305., 18246., 22008., 25379.,
29295., 33099., 36613., 40526., 43903., 47068., 50177.,
54459., 57496., 61353., 64993., 67824., 71757.],
[ 2182., 6286., 10273., 14212., 18125., 21717., 25602.,
29543., 32993., 36432., 39977., 43196., 47120., 50479.,
53572., 57702., 60687., 64071., 67382., 71237.],
[ 2065., 6144., 9957., 13940., 17881., 21754., 25439.,
29237., 32840., 36541., 39511., 43342., 46823., 50636.,
53757., 57248., 60289., 64381., 67124., 71203.],
[ 2040., 5963., 10038., 14009., 17905., 21711., 25475.,
29043., 32682., 36300., 39542., 42924., 46708., 50511.,
54318., 57356., 60479., 63775., 67051., 70676.]])
In [34]:
corr2d=(4.0*dd2d+rr2d*1.0-2.0*dr2d-2.0*rd2d)/rr2d
In [36]:
corr2d
Out[36]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [ ]:
In [2]:
correl200k0=np.array([ 4.28338704, 0.83319115, 0.5812655 , 0.45728027, 0.43167961,
0.39402051, 0.34946067, 0.33124462, 0.34585462, 0.33795129])
In [3]:
bins=np.array([ 0.002, 0.004, 0.006, 0.008, 0.01 , 0.012, 0.014, 0.016,
0.018, 0.02 ])
In [5]:
plt.plot(bins,correl200k0,'bo-')
Out[5]:
In [26]:
from __future__ import division
def xi0(dz,gamma,beta,kz):
return kz*(dz)**(-gamma)*(1.0+2.0*(1.0-gamma)/(3.0-gamma)*beta+(3.0-6.0*gamma+gamma*(2.0+gamma))/((3.0-gamma)*(5.0-gamma))*beta**2)
In [ ]:
In [ ]:
In [ ]:
In [23]:
xi0v=xi0(bins,1.0,2.0,1.0)
In [24]:
print xi0v
In [16]:
bins
Out[16]:
In [17]:
xi0v=xi0(0.002,1.0,2.0,1.0)
In [18]:
print xi0v
In [25]:
ans=curve_fit(xi0,bins,correl200k0)
In [26]:
ans
Out[26]:
In [33]:
perr = np.sqrt(np.diag(ans[1]))
In [34]:
perr
Out[34]:
In [38]:
plt.plot(bins, xi0(bins,1.74549968, 0.73620809, 0.09657359),'bo-')
plt.plot(bins,correl200k0,'r*-')
Out[38]:
In [39]:
ans1=curve_fit(xi0,bins,correl200k0,method='lm')
In [40]:
ans1
Out[40]:
In [ ]:
In [ ]:
In [ ]:
In [11]:
def musq(y,dz,zdth):
return 1.0/(1.0+y**2*(dz/zdth[:,None])**2)
In [12]:
def f(y,dz,zdth,beta,gamma):
musqv=musq(y,dz,zdth)
return 1.0+2.0*(1.0-gamma*musqv)/(3-gamma)*beta+(3.0-6.0*gamma*musqv+gamma*(2.0+gamma)*musqv**2)*beta**2/((3.0-gamma)*(5.0-gamma))
In [13]:
def xi(dzzdth,Kz,y,beta,gamma):
dz=dzzdth
zdth=dzzdth
musqv=musq(y,dz,zdth)
fv=f(y,dz,zdth,beta,gamma)
ans=Kz*musqv**(gamma/2.0)*fv/(dz)**gamma
return ans.ravel()
In [14]:
def xibetaf(dzzdth,Kz,y,gamma):
beta=1.0
return xi(dzzdth,Kz,y,beta,gamma)
In [17]:
ans=curve_fit(xibetaf,dzbin,corrells.ravel())
In [18]:
ans
Out[18]:
In [19]:
sol=ans[0]
In [22]:
def xi2d(dzzdth,Kz,y,beta,gamma):
dz=dzzdth
zdth=dzzdth
musqv=musq(y,dz,zdth)
fv=f(y,dz,zdth,beta,gamma)
ans=Kz*musqv**(gamma/2.0)*fv/(dz)**gamma
return ans
def xibetaf2d(dzzdth,Kz,y,gamma):
beta=1.0
return xi2d(dzzdth,Kz,y,beta,gamma)
In [23]:
plt.contour(xibetaf2d(dzbin,*sol))
Out[23]:
In [24]:
plt.contour(corrells)
Out[24]:
In [38]:
correlshr4=correlshr[0:4,0:4]
In [44]:
ans=curve_fit(xibetaf,dzbin[0:4],correlshr4.ravel())
In [45]:
ans
Out[45]:
In [ ]:
In [ ]:
In [39]:
dzbin[0:4]
Out[39]: