Fixed whitification. Numpy was not viewing the arrays but copying them
This commit is contained in:
parent
0aeb0d80c3
commit
ca379c91ee
@ -136,8 +136,12 @@ def whitify(density, L, cosmo, supergenerate=1, func='HU_WIGGLES'):
|
||||
if supergenerate > 1:
|
||||
|
||||
cond=np.isnan(density_hat_super)
|
||||
x = np.random.randn(np.count_nonzero(cond),2)/np.sqrt(2.0)
|
||||
density_hat_super[cond] = x[:,0] + 1j * x[:,1]
|
||||
print np.where(np.isnan(density_hat_super))[0].size
|
||||
Nz = np.count_nonzero(cond)
|
||||
density_hat_super.real[cond] = np.random.randn(Nz)
|
||||
density_hat_super.imag[cond] = np.random.randn(Nz)
|
||||
density_hat_super[cond] /= np.sqrt(2.0)
|
||||
print np.where(np.isnan(density_hat_super))[0].size
|
||||
|
||||
# Now we have to fix the Nyquist plane
|
||||
hNs = Ns/2
|
||||
@ -145,6 +149,8 @@ def whitify(density, L, cosmo, supergenerate=1, func='HU_WIGGLES'):
|
||||
Nplane = nyquist.size
|
||||
nyquist.flat[:Nplane/2] = np.sqrt(2.0)*nyquist.flat[Nplane:Nplane/2:-1].conj()
|
||||
|
||||
print np.where(np.isnan(density_hat_super))[0].size
|
||||
|
||||
return np.fft.irfftn(density_hat_super)*Ns**1.5
|
||||
|
||||
|
||||
|
@ -10,26 +10,34 @@ cosmo['omega_k_0'] = 0
|
||||
cosmo['omega_B_0']=0.049
|
||||
cosmo['SIGMA8']=0.8344
|
||||
cosmo['ns']=0.9624
|
||||
N0=128
|
||||
N0=256
|
||||
|
||||
doSimulation=True
|
||||
simShift=True
|
||||
|
||||
snap_id=int(sys.argv[1])
|
||||
astart=1/100.
|
||||
|
||||
if doSimulation:
|
||||
s = ct.loadRamsesAll("/nethome/lavaux/remote2/borgsim/", snap_id, doublePrecision=True)
|
||||
s = ct.loadRamsesAll("/nethome/lavaux/remote2/borgsim2/", snap_id, doublePrecision=True)
|
||||
astart=s.getTime()
|
||||
L = s.getBoxsize()
|
||||
|
||||
dsim = ct.cicParticles(s.getPositions(), L, N0)
|
||||
p = s.getPositions()
|
||||
Nsim = int( np.round( p[0].size**(1./3)) )
|
||||
print("Nsim = %d" % Nsim)
|
||||
|
||||
if simShift:
|
||||
p = [(q-0.5*L/Nsim)%L for q in p]
|
||||
|
||||
dsim = ct.cicParticles(p[::-1], L, N0)
|
||||
dsim /= np.average(np.average(np.average(dsim, axis=0), axis=0), axis=0)
|
||||
dsim -= 1
|
||||
|
||||
dsim_hat = np.fft.rfftn(dsim)*(L/N0)**3
|
||||
Psim, bsim = bic.bin_power(np.abs(dsim_hat)**2/L**3, L, range=(0,1.), bins=150)
|
||||
|
||||
pos,_,density,N,L,_,_ = bic.run_generation("initial_density_1380.dat", 0.001, astart, cosmo, supersample=2, do_lpt2=True)
|
||||
pos,_,density,N,L,_,_ = bic.run_generation("initial_density_2588.dat", 0.001, astart, cosmo, supersample=2, do_lpt2=True)
|
||||
|
||||
dcic = ct.cicParticles(pos, L, N0)
|
||||
dcic /= np.average(np.average(np.average(dcic, axis=0), axis=0), axis=0)
|
||||
@ -50,7 +58,7 @@ Pref, bref = bic.compute_ref_power(L, N0, cosmo, range=(0,1.), bins=150)
|
||||
Pcic /= D1_0**2
|
||||
Pdens /= D1_0**2
|
||||
|
||||
borg_evolved = ct.read_borg_vol("final_density_1380.dat")
|
||||
dborg_hat = np.fft.rfftn(borg_evolved.density)*L**3/borg_evolved.density.size
|
||||
#borg_evolved = ct.read_borg_vol("final_density_1380.dat")
|
||||
#dborg_hat = np.fft.rfftn(borg_evolved.density)*L**3/borg_evolved.density.size
|
||||
|
||||
Pborg, bborg = bic.bin_power(np.abs(dborg_hat)**2/L**3, L, range=(0,1.),bins=150)
|
||||
#Pborg, bborg = bic.bin_power(np.abs(dborg_hat)**2/L**3, L, range=(0,1.),bins=150)
|
||||
|
Loading…
Reference in New Issue
Block a user