<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 Import NumPy and SciPy (not needed when using --pylab)
Step1:  Load data from file
Step2:  Construct histogram from data
There are several histrogram commands
Step3:  Interpolate histogram output -> p(z); n.b. that you can also use numerical quadrature to get $P(z)$ directly.
Step4:  Use numerical integration to get $P(z) = \int_0^\infty p(z') dz'$
<ASSISTANT_TASK:>
Python Code:
%pylab inline
import scipy as sp
Explanation: Import NumPy and SciPy (not needed when using --pylab)
End of explanation
zz = np.loadtxt('wiggleZ_DR1_z.dat',dtype='float'); # Load WiggleZ redshifts
np.min(zz) # Check bounds
np.max(zz)
Explanation: Load data from file
End of explanation
nbins = 50; # Is this a good choice? 
n, bins, patches = hist(zz,nbins)
x = bins[0:nbins] + (bins[2]-bins[1])/2; # Convert bin edges to centres, chopping the last
Explanation: Construct histogram from data
There are several histrogram commands: hist() will be fine here, but note the syntax below. Also note that the bin edges are returned, so that there will nbins+1 of these.
End of explanation
# Import the function you need
from scipy.interpolate import interp1d
# Build an interpolation function for p(z) that accepts an arbitrary redshift z
zdist = interp1d(x,n, kind='cubic', bounds_error = False, fill_value=0)
z = linspace(0,2,100); plot(z,zdist(z)) # Test your interpolation function out
Explanation: Interpolate histogram output -> p(z); n.b. that you can also use numerical quadrature to get $P(z)$ directly.
End of explanation
# Import the function you need
from scipy import integrate
intPz = lambda z0: integrate.quad(zdist,0,z0)[0]
total = intPz(5) # Get normalisation constant by evaluating P(z->\infty)
total # Check that this worked
Explanation: Use numerical integration to get $P(z) = \int_0^\infty p(z') dz'$
End of explanation
<END_TASK>