Merge branch 'master' of bitbucket.org:glavaux/cosmotool
This commit is contained in:
commit
adf14da4b4
@ -7,8 +7,9 @@ include(GetGitRevisionDescription)
|
||||
include(ExternalProject)
|
||||
include(FindOpenMP)
|
||||
|
||||
include(FindPythonLibs)
|
||||
include(FindPythonInterp)
|
||||
set(NumPy_FIND_REQUIRED TRUE)
|
||||
include(${CMAKE_SOURCE_DIR}/FindPyLibs.cmake)
|
||||
include(${CMAKE_SOURCE_DIR}/FindNumPy.cmake)
|
||||
|
||||
|
||||
|
30
FindPyLibs.cmake
Normal file
30
FindPyLibs.cmake
Normal file
@ -0,0 +1,30 @@
|
||||
execute_process(COMMAND "${PYTHON_EXECUTABLE}" "-c"
|
||||
"import distutils.sysconfig as cs; import os; v=cs.get_config_vars(); print(os.path.join(v['LIBDIR'],v['LDLIBRARY']));"
|
||||
RESULT_VARIABLE _PYLIB_SEARCH_SUCCESS
|
||||
OUTPUT_VARIABLE _PYLIB_VALUES_OUTPUT
|
||||
ERROR_VARIABLE _PYLIB_ERROR_VALUE
|
||||
OUTPUT_STRIP_TRAILING_WHITESPACE)
|
||||
|
||||
execute_process(COMMAND "${PYTHON_EXECUTABLE}" "-c"
|
||||
"import distutils.sysconfig as cs; import os; v=cs.get_config_vars(); print(v['INCLUDEPY']);"
|
||||
RESULT_VARIABLE _PYINC_SEARCH_SUCCESS
|
||||
OUTPUT_VARIABLE _PYINC_VALUES_OUTPUT
|
||||
ERROR_VARIABLE _PYINC_ERROR_VALUE
|
||||
OUTPUT_STRIP_TRAILING_WHITESPACE)
|
||||
|
||||
|
||||
if(NOT _PYLIB_SEARCH_SUCCESS MATCHES 0)
|
||||
message(FATAL_ERROR
|
||||
"PyLib search failure:\n${_PYLIB_ERROR_VALUE}")
|
||||
return()
|
||||
endif()
|
||||
|
||||
if(NOT _PYINC_SEARCH_SUCCESS MATCHES 0)
|
||||
message(FATAL_ERROR
|
||||
"PyInc search failure:\n${_PYINC_ERROR_VALUE}")
|
||||
return()
|
||||
endif()
|
||||
|
||||
|
||||
set(PYTHON_LIBRARY ${_PYLIB_VALUES_OUTPUT} CACHE PATH "Python runtime library path")
|
||||
set(PYTHON_INCLUDE_PATH ${_PYINC_VALUES_OUTPUT} CACHE PATH "Python runtime include path")
|
@ -1 +1,2 @@
|
||||
from _cosmotool import *
|
||||
from borg import read_borg_vol
|
||||
|
265
python/cosmotool/borg.py
Normal file
265
python/cosmotool/borg.py
Normal file
@ -0,0 +1,265 @@
|
||||
###
|
||||
### BORG code is from J. Jasche
|
||||
###
|
||||
import StringIO
|
||||
import numpy as np
|
||||
from numpy import *
|
||||
import os.path
|
||||
import array
|
||||
import glob
|
||||
|
||||
class BorgVolume(object):
|
||||
|
||||
def __init__(self, density, ranges):
|
||||
self._density = density
|
||||
self._ranges = ranges
|
||||
|
||||
property density:
|
||||
def __get__(self):
|
||||
return self._density
|
||||
|
||||
property ranges:
|
||||
def __get__(self):
|
||||
return self._ranges
|
||||
|
||||
def build_filelist(fdir):
|
||||
#builds list of all borg density fields which may be distributed over several directories
|
||||
|
||||
fname_0=glob.glob(fdir[0]+'initial_density_*')
|
||||
fname_1=glob.glob(fdir[0]+'final_density_*')
|
||||
|
||||
fdir=fdir[1:] #eliminate first element
|
||||
|
||||
for fd in fdir:
|
||||
fname_0=fname_0+glob.glob(fd+'initial_density_*')
|
||||
fname_1=fname_1+glob.glob(fd+'final_density_*')
|
||||
|
||||
|
||||
return fname_0, fname_1
|
||||
|
||||
def read_borg_vol(BORGFILE):
|
||||
""" Reading routine for BORG data
|
||||
"""
|
||||
|
||||
openfile=open(BORGFILE,'rb')
|
||||
|
||||
period=0
|
||||
N0=0
|
||||
N1=0
|
||||
N2=0
|
||||
|
||||
xmin=0
|
||||
xmax=0
|
||||
|
||||
ymin=0
|
||||
ymax=0
|
||||
|
||||
zmin=0
|
||||
zmax=0
|
||||
nlines=0
|
||||
|
||||
while True:
|
||||
line=openfile.readline()
|
||||
s=line.rstrip('\n')
|
||||
r=s.rsplit(' ')
|
||||
|
||||
if size(r)==5 :
|
||||
if r[0] =="define":
|
||||
if r[1]=="Lattice" : N0=int(r[2])
|
||||
N1=int(r[3])
|
||||
N2=int(r[4])
|
||||
|
||||
|
||||
if size(r)==11 :
|
||||
if r[4] =="BoundingBox": xmin=float(r[5])
|
||||
xmax=float(r[6])
|
||||
ymin=float(r[7])
|
||||
ymax=float(r[8])
|
||||
zmin=float(r[9])
|
||||
zmax=float(r[10].rstrip(','))
|
||||
|
||||
if r[0]=='@1': break
|
||||
|
||||
ranges=[]
|
||||
ranges.append(xmin)
|
||||
ranges.append(xmax)
|
||||
ranges.append(ymin)
|
||||
ranges.append(ymax)
|
||||
ranges.append(zmin)
|
||||
ranges.append(zmax)
|
||||
|
||||
#now read data
|
||||
data=np.fromfile(openfile, '>f4')
|
||||
data=data.reshape(N2,N0,N1)
|
||||
return BorgVolume(data,ranges)
|
||||
|
||||
def read_spec( fname ):
|
||||
""" Reading routine for ARES spectrum samples
|
||||
"""
|
||||
x,y=np.loadtxt( fname ,usecols=(0,1),unpack=True)
|
||||
|
||||
return x , y
|
||||
|
||||
|
||||
def read_bias_nmean( fname ):
|
||||
""" Reading routine for ARES bias data
|
||||
"""
|
||||
x,b0,b1,nmean=np.loadtxt( fname ,usecols=(0,1,2,3),unpack=True)
|
||||
|
||||
return x , b0, b1, nmean
|
||||
|
||||
def read_nmean( fname ):
|
||||
""" Reading routine for BORG bias data
|
||||
"""
|
||||
x,nmean=np.loadtxt( fname ,usecols=(0,1),unpack=True)
|
||||
|
||||
return x, nmean
|
||||
|
||||
def get_grid_values(xx,data, ranges):
|
||||
""" return values at grid positions
|
||||
"""
|
||||
xmin=ranges[0]
|
||||
xmax=ranges[1]
|
||||
|
||||
ymin=ranges[2]
|
||||
ymax=ranges[3]
|
||||
|
||||
zmin=ranges[4]
|
||||
zmax=ranges[5]
|
||||
|
||||
Lx= xmax-xmin
|
||||
Ly= ymax-ymin
|
||||
Lz= zmax-zmin
|
||||
|
||||
Nx=shape(data)[0]
|
||||
Ny=shape(data)[1]
|
||||
Nz=shape(data)[2]
|
||||
|
||||
dx=Lx/float(Nx)
|
||||
dy=Ly/float(Ny)
|
||||
dz=Lz/float(Nz)
|
||||
|
||||
idx=(xx[:,0]-xmin)/dx
|
||||
idy=(xx[:,1]-ymin)/dz
|
||||
idz=(xx[:,2]-zmin)/dy
|
||||
|
||||
idx=idx.astype(int)
|
||||
idy=idy.astype(int)
|
||||
idz=idz.astype(int)
|
||||
|
||||
idflag=np.where( (idx>-1)*(idx<Nx)*(idy>-1)*(idy<Ny)*(idz>-1)*(idz<Nz) )
|
||||
|
||||
flag=[False]*len(xx)
|
||||
vals=[-999.]*len(xx)
|
||||
|
||||
flag=np.array(flag)
|
||||
vals=np.array(vals)
|
||||
flag[idflag]=True
|
||||
vals[idflag]=data[idx[idflag],idy[idflag],idz[idflag]]
|
||||
|
||||
return vals,flag
|
||||
|
||||
def get_mean_density(fdir, smin, step):
|
||||
""" estimate ensemble mean
|
||||
"""
|
||||
print '-'*60
|
||||
print 'Get 3D ensemble mean density field'
|
||||
print '-'*60
|
||||
|
||||
fname0 = fdir + 'initial_density_'+str(0)+'.dat'
|
||||
fname1 = fdir + 'final_density_'+str(0)+'.dat'
|
||||
|
||||
MEAN0,ranges=read_borg_vol(fname0)
|
||||
MEAN0=MEAN0*0.;
|
||||
VAR0=copy(MEAN0)
|
||||
MEAN1=copy(MEAN0)
|
||||
VAR1=copy(MEAN0)
|
||||
norm=0.
|
||||
|
||||
idat=smin
|
||||
|
||||
fname0 = fdir + 'initial_density_'+str(idat)+'.dat'
|
||||
fname1 = fdir + 'final_density_'+str(idat)+'.dat'
|
||||
#and (idat<smin+1000)
|
||||
while((os.path.exists(fname0))):
|
||||
auxdata0,auxranges0=read_borg_vol(fname0)
|
||||
auxdata1,auxranges1=read_borg_vol(fname1)
|
||||
auxx0=auxdata0
|
||||
auxx1=auxdata1
|
||||
MEAN0+=auxx0
|
||||
VAR0+=auxx0**2
|
||||
MEAN1+=auxx1
|
||||
VAR1+=auxx1**2
|
||||
|
||||
norm+=1
|
||||
idat+=step
|
||||
fname0 = fdir + 'initial_density_'+str(idat)+'.dat'
|
||||
fname1 = fdir + 'final_density_'+str(idat)+'.dat'
|
||||
del auxranges0
|
||||
del auxdata0
|
||||
del auxranges1
|
||||
del auxdata1
|
||||
|
||||
MEAN0/=norm
|
||||
VAR0/=norm
|
||||
VAR0-=MEAN0**2
|
||||
VAR0=sqrt(fabs(VAR0))
|
||||
|
||||
MEAN1/=norm
|
||||
VAR1/=norm
|
||||
VAR1-=MEAN1**2
|
||||
VAR1=sqrt(fabs(VAR1))
|
||||
|
||||
|
||||
return MEAN0,VAR0,MEAN1,VAR1,ranges
|
||||
|
||||
def get_mean_density_fdir(fdir,init,steps):
|
||||
""" estimate ensemble mean
|
||||
"""
|
||||
print '-'*60
|
||||
print 'Get 3D ensemble mean density field'
|
||||
print '-'*60
|
||||
|
||||
fname0,fname1=build_filelist(fdir)
|
||||
|
||||
fname0=fname0[init::steps]
|
||||
fname1=fname1[init::steps]
|
||||
|
||||
MEAN0,ranges=read_borg_vol(fname0[0])
|
||||
MEAN0=MEAN0*0.;
|
||||
VAR0=copy(MEAN0)
|
||||
MEAN1=copy(MEAN0)
|
||||
VAR1=copy(MEAN0)
|
||||
norm0=0.
|
||||
norm1=0.
|
||||
|
||||
for fn in fname0:
|
||||
auxdata0,auxranges0=read_borg_vol(fn)
|
||||
MEAN0+=auxdata0
|
||||
VAR0+=auxdata0**2.
|
||||
norm0+=1.
|
||||
del auxranges0
|
||||
del auxdata0
|
||||
|
||||
for fn in fname1:
|
||||
auxdata1,auxranges1=read_borg_vol(fn)
|
||||
MEAN1+=auxdata1
|
||||
VAR1+=auxdata1**2.
|
||||
norm1+=1.
|
||||
del auxranges1
|
||||
del auxdata1
|
||||
|
||||
MEAN0/=norm0
|
||||
VAR0/=norm0
|
||||
VAR0-=MEAN0**2
|
||||
VAR0=sqrt(fabs(VAR0))
|
||||
|
||||
MEAN1/=norm1
|
||||
VAR1/=norm1
|
||||
VAR1-=MEAN1**2
|
||||
VAR1=sqrt(fabs(VAR1))
|
||||
|
||||
|
||||
return MEAN0,VAR0,MEAN1,VAR1,ranges
|
||||
|
||||
|
Loading…
Reference in New Issue
Block a user