Modifying 3D Correlation Function

I'm changing the 3D correlation function code to accept data and randoms (similar to what I just did for the 2D correlation function code).

This is proving to be a little tricky because in the current version I do some tricks to make sure the periodic boundary conditions don't add power incorrectly, and I also re-scale the data because the correlation function is expecting it to be in a 1X1X1 box.

After talking with Alexia, I've decided to do these conversions outside of the correlation function in python, such that the code expects the data to be in a 1x1x1 box.

Results to come soon.....

Matching Results!

It turned out I had was printing out the same random catalog twice from the old code! So now everything matches!

Old:
# theta (degrees) crossw(theta)
1.1250000000 9.9488993345
1.3750000000 2.3509358289
1.6250000000 1.6241751936
1.8750000000 1.2341323430
2.1250000000 1.1220636000
2.3750000000 0.9348225839
2.6250000000 0.8468507873
2.8750000000 0.8237087312
3.1250000000 0.8176693767
3.3750000000 0.7284754936
3.6250000000 0.7714957667
3.8750000000 0.7788989999
4.1250000000 0.8089800491
4.3750000000 0.7564553032
4.6250000000 0.6803548995
4.8750000000 0.6352504229
5.1250000000 0.6443996748
5.3750000000 0.5349427778
5.6250000000 0.5059338935
5.8750000000 0.4477430500

New:
# theta (degrees) crossw(theta)
1.1250000000 9.9488993345
1.3750000000 2.3509358289
1.6250000000 1.6241751936
1.8750000000 1.2341323430
2.1250000000 1.1220636000
2.3750000000 0.9348225839
2.6250000000 0.8468507873
2.8750000000 0.8237087312
3.1250000000 0.8176693767
3.3750000000 0.7284754936
3.6250000000 0.7714957667
3.8750000000 0.7788989999
4.1250000000 0.8089800491
4.3750000000 0.7564553032
4.6250000000 0.6803548995
4.8750000000 0.6352504229
5.1250000000 0.6443996748
5.3750000000 0.5349427778
5.6250000000 0.5059338935
5.8750000000 0.4477430500

WPS Details Compairison

Old Code:
bin dd dr rd rr
0 43185.000000 2336.000000 2352.600000 5044.200000
1 4494.000000 1037.400000 1059.200000 1966.600000
2 3790.000000 1210.200000 1214.800000 2273.400000
3 3392.000000 1400.200000 1389.800000 2550.000000
4 3556.000000 1594.000000 1681.000000 2794.000000
5 3506.000000 1823.200000 1922.600000 3108.000000
6 3564.000000 2001.400000 2054.600000 3351.400000
7 3852.000000 2182.000000 2314.000000 3557.200000
8 4056.000000 2389.400000 2425.200000 3827.600000
9 4050.000000 2566.000000 2612.400000 4035.600000
10 4614.000000 2760.600000 2791.800000 4375.000000
11 4968.000000 2998.600000 3088.800000 4575.200000
12 5392.000000 3139.600000 3286.800000 4800.400000
13 5542.000000 3347.000000 3566.800000 5038.600000
14 5544.000000 3551.600000 3687.600000 5257.200000
15 5600.000000 3746.600000 3873.400000 5582.400000
16 5946.000000 3924.800000 4025.800000 5795.600000
17 5452.000000 4065.800000 4100.200000 6035.400000
18 5464.000000 4235.600000 4345.200000 6207.400000
19 5412.000000 4380.200000 4527.400000 6357.600000

New:
0 43185.000000 2352.600000 2352.600000 5442.800000
1 4494.000000 1059.200000 1059.200000 2048.600000
2 3790.000000 1214.800000 1214.800000 2347.800000
3 3392.000000 1389.800000 1389.800000 2611.800000
4 3556.000000 1681.000000 1681.000000 2922.400000
5 3506.000000 1922.600000 1922.600000 3123.200000
6 3564.000000 2054.600000 2054.600000 3369.000000
7 3852.000000 2314.000000 2314.000000 3661.200000
8 4056.000000 2425.200000 2425.200000 3896.600000
9 4050.000000 2612.400000 2612.400000 4160.200000
10 4614.000000 2791.800000 2791.800000 4406.000000
11 4968.000000 3088.800000 3088.800000 4704.600000
12 5392.000000 3286.800000 3286.800000 4951.800000
13 5542.000000 3566.800000 3566.800000 5215.400000
14 5544.000000 3687.600000 3687.600000 5396.000000
15 5600.000000 3873.400000 3873.400000 5598.200000
16 5946.000000 4025.800000 4025.800000 5856.000000
17 5452.000000 4100.200000 4100.200000 5987.400000
18 5464.000000 4345.200000 4345.200000 6225.600000
19 5412.000000 4527.400000 4527.400000 6421.800000

Looks like DD and RD match between old and new.... and in the new the DR and RD are the same, so I must be pointing to the same random file twice in the new code, better check that!

Comparing Code

I've got the old correlation function code (where it creates the randoms inside the code) spitting out the randoms it is using to a file, and I am then comparing the correlation function to my new code using these random files. Here are the two correlation function values:

New code:
1.1250000000 8.0698537517
1.3750000000 2.1596212047
1.6250000000 1.5794360678
1.8750000000 1.2344743089
2.1250000000 1.0663837941
2.3750000000 0.8913934426
2.6250000000 0.8381715643
2.8750000000 0.7880476347
3.1250000000 0.7961299595
3.3750000000 0.7176097303
3.6250000000 0.7799364503
3.8750000000 0.7428899375
4.1250000000 0.7613797003
4.3750000000 0.6948268589
4.6250000000 0.6606375093
4.8750000000 0.6165195956
5.1250000000 0.6404371585
5.3750000000 0.5409693690
5.6250000000 0.4817527628
5.8750000000 0.4327447133

Old Code:
1.1250000000 8.6318147575
1.3750000000 2.2190582732
1.6250000000 1.6004222750
1.8750000000 1.2360784314
2.1250000000 1.1005726557
2.3750000000 0.9228442728
2.6250000000 0.8531956794
2.8750000000 0.8189587316
3.1250000000 0.8018079214
3.3750000000 0.7203885420
3.6250000000 0.7855085714
3.8750000000 0.7553331002
4.1250000000 0.7845179568
4.3750000000 0.7277418330
4.6250000000 0.6775469832
4.8750000000 0.6381484666
5.1250000000 0.6541169163
5.3750000000 0.5503197800
5.6250000000 0.4978896156
5.8750000000 0.4501698754

The numbers are close, but not exactly the same... they should be exactly the same I think. :( I'm looking deeper into the calculation to see where the difference is exactly.

Masking both Data and Randoms

I've masked both the data and the randoms now:


The randoms (magenta) now fall exactly on top of the data (cyan). The code to do this is as follows:

;Read in data
file1 = "/home/shirleyho/research/SDSS_PK/powerspec_dir/DATA/ra_dec_z_run_rerun_camcol.z_gt_0.01"

readcol, file1, RA, DEC, z, RUN, rerun, camcol, x, x, x, x, x, x, x, x, format='(F,F,F,D,D,D,D,D,F,F,F,F,F,F)'

datasize = n_elements(ra)

;Read in mask
file = "/home/shirleyho/research/SDSS_PK/powerspec_dir/MASK_SDSS/dir_fits/tmp_unified_score_ebv_psf_bsm_nside1024.fits"

read_fits_map, file, weight, nside=nside_in, ordering=order


index = 0L
;Convert ra/dec to phi/theta expected in HEALPix
dataphi = ra*!pi/180.
datatheta = (90.0-dec)*!pi/180.
;Will be filled in with weight corresponding to ra/dec of data
dataweight = findgen(datasize)*0.0

;Function I wrote to get the weight for the data
junk = maskdata(datatheta, dataphi, index, datasize, nside_in, weight, dataweight)
;Only keeps data with a weight greater than 0.2 (inside mask)
datacut = where(dataweight GE 0.2)

;Make the randoms
randomsize = datasize*5L
randomtheta = findgen(randomsize)*0.0
randomphi = findgen(randomsize)*0.0
index = 0L
;Function I wrote to fill in randomtheta and randomphi with points inside the mask
junk = makemaskrandoms(randomtheta, randomphi, index, randomsize, nside_in, weight)

;Plot the data and the randoms
plot, ra[sub], dec[sub], psym=3, symsize=2, xr=[0,360],yr=[0,90], XTITLE = "RA", YTITLE = "DEC", TITLE = "RA vs Dec of Data and Randoms", charsize = 2, charthick = 1, thick = 2, xthick=2, ythick=2
oplot, ra[datacut], dec[datacut], psym=3, color=fsc_color('cyan')
oplot, randomra[0:100000], randomdec[0:100000], psy = 3, color=fsc_color('magenta')

Now to compute some correlation functions!

Plot of Data / Randoms

Here is a plot of the data (cyan) and randoms (magenta) using the masks, and the method outlined in this post:

I'm a little confused as to why some of the data isn't represented in the randoms. It seems the mask is different than the data, perhaps this is because there is less completeness in the regions that aren't covered in the mask. I'll ask Shirley about this, but in the mean time, I'll apply the mask to both the data and the randoms, to insure the area covered are the same.


11/2 Update: Heard back from Shirley she says:

"This is deliberate. And your plot reflects it :)

There is a reason to cut out certain areas in the mask, and I never said that the data covers the same area as the mask :) ! I am testing things right now for different areas to cut out :) ! (due to extinction, stars and possible photometric calibration problems)"

So I looks like I should just apply the mask to both the data and the randoms before doing the cross correlation, as I suspected.