Tuesday, May 10, 2011

Different Reconstruction

Here is the reconstruction I ran last night with a different photometric distribution function, below shows the photometric and spectroscopic distributions:

I used the same number of bins and angles as yesterdays post. Here are the correlation functions:

And here are the reconstructions with several different binnings. None of them match as well as yesterdays post, but both 40 and 39 are close-ish:

I talked to Alexia this morning, she suggested that I do the binning sensitivity check and also look at how xmpl_rv (in getximat) compares to sximat. I'm doing that now. Most posts to come shortly.

Here's the code to make the following (also in logs/110509log.py):
import numpy as N
from pylab import *
from matplotlib import *
import crlold as crl
from correlationData import *
import galsphere as gals
import datetime as dt

workingDir = 'run201159_123'

runConstantsFile = './'+workingDir+'/runConstants.dat'


rbinspace = int((rbinmax-rbinmin)/(nrbins))
rbins = makeRbins(rbinmin,rbinmax,rbinspace)

#read in masks

#Import the photometric redshift distribution
photoCD = readPhotoCDfile(photoCDfile)

#Create a matrix to input the correlation function data
# Write wps Matrix to file
#The number of redshift bins you want to skip for the
#reconstruction (this removes noisy bins from calculation)
skip = 0
makeWPSmatrix2(workingDir,rbins,corrBins, wpsMatrixFile, skip)

#Create a matrix to input the correlation function data
# Write Xi Matrix to file
makeXimatrix2(workingDir,rbins,nxicorr,xiMatrixFile, skip)

#Plot wps correlation functions
numtoplot = int(nrbins-skip) #number of correlation functions you want to display on plot
plotWPS(skip, numtoplot, workingDir, rbins)
numtoplot = 10
plotWPSScreen(skip, numtoplot, workingDir, rbins)

#Plot xi correlation functions
plotXi(skip, numtoplot, workingDir, rbins)
plotXiScreen(skip, numtoplot, workingDir, rbins)

#from commandsAlexia import *
xissfile = xiMatrixFile
wpsfile = wpsMatrixFile
nphibins = size(rbins)

# Read in the data
#spectroscopic autocorrelation function
rs=X[0,:]/1000. #want box units (Gpc, not Mpc)
xispec=X[1:,:] #xi in nspecbins

#2D cross correlation between spec and photo samples
#rlos=getcomovingD(rlos) #adjusts from redshift to Mpc/h
rlos=rlos/1000 #want box units (Gpc, not Mpc)

#covariance matrix of wps measurements, wps vector arranged
covar = N.identity(wps.shape[0],float)
xierr = N.identity(xispec.shape[1],float)

# compute the xi matrix using the spec sample,
npbins=39#the number of bins you want reconstructed
reconstructmax = 0.5 #in Gpc
reconstructionmin = 0.0 #in Gpc
phibsize = (reconstructmax - reconstructionmin)/npbins

phibins = linspace(start=reconstructionmin,stop=reconstructmax,num=npbins,endpoint=True)

#phibins =(linspace(2,npbins,num=npbins,endpoint=False)+0.5)*(1./npbins)

# invert to recover phi
for i in range(len(phircvd.x)):
if (phircvd.x[i] < 0):

print 'the recovered selection function is not normalized.'

outphi4 = phircvd.x
phibins4 = phibins

start = 1
end = npbins - 1

#junkbins = linspace(start=0,stop=0.5,num=100,endpoint=True)
junkbins = phibins
junk = histogram(photoCD/1000., bins=junkbins)
junkdata = junk[0]
plot(junkbins[start:end], junkdata[start:end], 'r.-', label = 'photometric data')
plot((phibins4[start:end]+1.0*phibsize), outphi4[start:end]*sum(junkdata[start:end])/sum(outphi4[start:end]), 'g-', label='reconstruction %i'%npbins)

xlabel('comoving distance (Gpc/h)')
ylabel('number of objects per bin')
title('Redshift Distribution')
from pylab import *

No comments:

Post a Comment