Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add compute method and remove Catalog interface from CatalogMesh. #498

Merged
merged 4 commits into from
Jun 1, 2018
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
62 changes: 31 additions & 31 deletions nbodykit/algorithms/convpower.py
Original file line number Diff line number Diff line change
Expand Up @@ -144,9 +144,9 @@ def __init__(self, first, poles,
use_fkp_weights=False,
P0_FKP=None):

first = _cast_source(first, Nmesh=Nmesh)
first = _cast_mesh(first, Nmesh=Nmesh)
if second is not None:
second = _cast_source(second, Nmesh=Nmesh)
second = _cast_mesh(second, Nmesh=Nmesh)
else:
second = first

Expand Down Expand Up @@ -194,22 +194,22 @@ def __init__(self, first, poles,

# add FKP weights
if use_fkp_weights:
for source in [self.first, self.second]:
for mesh in [self.first, self.second]:
if self.comm.rank == 0:
args = (source.fkp_weight, P0_FKP)
args = (mesh.fkp_weight, P0_FKP)
self.logger.info("adding FKP weights as the '%s' column, using P0 = %.4e" %args)

for name in ['data', 'randoms']:
for name in ['data', 'randoms']:

# print a warning if we are overwriting a non-default column
old_fkp_weights = source[name][source.fkp_weight]
if source.compute(old_fkp_weights.sum()) != len(old_fkp_weights):
warn = "it appears that we are overwriting FKP weights for the '%s' " %name
warn += "source in FKPCatalog when using 'use_fkp_weights=True' in ConvolvedFFTPower"
warnings.warn(warn)
# print a warning if we are overwriting a non-default column
old_fkp_weights = mesh.source[name][mesh.fkp_weight]
if mesh.source.compute(old_fkp_weights.sum()) != len(old_fkp_weights):
warn = "it appears that we are overwriting FKP weights for the '%s' " %name
warn += "source in FKPCatalog when using 'use_fkp_weights=True' in ConvolvedFFTPower"
warnings.warn(warn)

nbar = source[name][source.nbar]
source[name][source.fkp_weight] = 1.0 / (1. + P0_FKP * nbar)
nbar = mesh.source[name][mesh.nbar]
mesh.source[name][mesh.fkp_weight] = 1.0 / (1. + P0_FKP * nbar)

# store meta-data
self.attrs = {}
Expand Down Expand Up @@ -461,7 +461,7 @@ class uses the spherical harmonic addition theorem such that
Ylms = [[get_real_Ylm(l,m) for m in range(-l, l+1)] for l in poles[1:]]

# paint the 1st FKP density field to the mesh (paints: data - alpha*randoms, essentially)
rfield1 = self.first.paint(Nmesh=self.attrs['Nmesh'])
rfield1 = self.first.compute(Nmesh=self.attrs['Nmesh'])
meta1 = rfield1.attrs.copy()
if rank == 0:
self.logger.info("%s painting of 'first' done" %self.first.window)
Expand All @@ -485,7 +485,7 @@ class uses the spherical harmonic addition theorem such that
if self.first is not self.second:

# paint the second field
rfield2 = self.second.paint(Nmesh=self.attrs['Nmesh'])
rfield2 = self.second.compute(Nmesh=self.attrs['Nmesh'])
meta2 = rfield2.attrs.copy()
if rank == 0: self.logger.info("%s painting of 'second' done" %self.second.window)

Expand Down Expand Up @@ -670,11 +670,11 @@ def normalization(self, name, alpha):
if name+'.norm' not in self.attrs:

# the selection (same for first/second)
sel = self.first.compute(self.first[name][self.first.selection])
sel = self.first.source.compute(self.first.source[name][self.first.selection])

# selected first/second meshes for "name" (data or randoms)
first = self.first[name][sel]
second = self.second[name][sel]
first = self.first.source[name][sel]
second = self.second.source[name][sel]

# these are assumed the same b/w first and second meshes
comp_weight = first[self.first.comp_weight]
Expand All @@ -690,7 +690,7 @@ def normalization(self, name, alpha):
A = nbar*comp_weight*fkp_weight1*fkp_weight2
if name == 'randoms':
A *= alpha
A = self.first.compute(A.sum())
A = first.compute(A.sum())
self.attrs[name+'.norm'] = self.comm.allreduce(A)

return self.attrs[name+'.norm']
Expand Down Expand Up @@ -718,11 +718,11 @@ def shotnoise(self, alpha):
for name in ['data', 'randoms']:

# the selection (same for first/second)
sel = self.first.compute(self.first[name][self.first.selection])
sel = self.first.source.compute(self.first.source[name][self.first.selection])

# selected first/second meshes for "name" (data or randoms)
first = self.first[name][sel]
second = self.second[name][sel]
first = self.first.source[name][sel]
second = self.second.source[name][sel]

# completeness weights (assumed same for first/second)
comp_weight = first[self.first.comp_weight]
Expand All @@ -745,26 +745,26 @@ def shotnoise(self, alpha):
# divide by normalization from randoms
return Pshot / self.attrs['randoms.norm']

def _cast_source(source, Nmesh):
def _cast_mesh(mesh, Nmesh):
"""
Cast an object to a MeshSource. Nmesh is used only on FKPCatalog
"""
from nbodykit.source.catalog import FKPCatalog
from nbodykit.source.catalogmesh import FKPCatalogMesh
if not isinstance(source, (FKPCatalogMesh, FKPCatalog)):
if not isinstance(mesh, (FKPCatalogMesh, FKPCatalog)):
raise TypeError("input sources should be a FKPCatalog or FKPCatalogMesh")

if isinstance(source, FKPCatalog):
if isinstance(mesh, FKPCatalog):
# if input is CatalogSource, use defaults to make it into a mesh
if not isinstance(source, FKPCatalogMesh):
source = source.to_mesh(Nmesh=Nmesh, dtype='f8', compensated=False)
if not isinstance(mesh, FKPCatalogMesh):
mesh = mesh.to_mesh(Nmesh=Nmesh, dtype='f8', compensated=False)

if Nmesh is not None and any(source.attrs['Nmesh'] != Nmesh):
raise ValueError(("Mismatched Nmesh between __init__ and source.attrs; "
if Nmesh is not None and any(mesh.attrs['Nmesh'] != Nmesh):
raise ValueError(("Mismatched Nmesh between __init__ and mesh.attrs; "
"if trying to re-sample with a different mesh, specify "
"`Nmesh` as keyword of to_mesh()"))

return source
return mesh

def get_compensation(mesh):
toret = None
Expand All @@ -784,7 +784,7 @@ def copy_meta(attrs, meta, prefix=""):

def is_valid_crosscorr(first, second):

if second.base is not first.base:
if second.source is not first.source:
return False

same_cols = ['selection', 'comp_weight', 'nbar']
Expand Down
12 changes: 6 additions & 6 deletions nbodykit/algorithms/fftpower.py
Original file line number Diff line number Diff line change
Expand Up @@ -97,13 +97,13 @@ def _compute_3d_power(self):
p3d : array_like (complex)
the 3D complex array holding the power spectrum
"""
c1 = self.first.paint(mode='complex', Nmesh=self.attrs['Nmesh'])
c1 = self.first.compute(mode='complex', Nmesh=self.attrs['Nmesh'])

# compute the auto power of single supplied field
if self.first is self.second:
c2 = c1
else:
c2 = self.second.paint(mode='complex', Nmesh=self.attrs['Nmesh'])
c2 = self.second.compute(mode='complex', Nmesh=self.attrs['Nmesh'])

# calculate the 3d power spectrum, slab-by-slab to save memory
p3d = c1
Expand Down Expand Up @@ -356,13 +356,13 @@ def _compute_3d_power(self):
p3d : array_like (complex)
the 3D complex array holding the power spectrum
"""
c1 = self.first.paint(mode='complex', Nmesh=self.attrs['Nmesh'])
c1 = self.first.compute(mode='complex', Nmesh=self.attrs['Nmesh'])

# compute the auto power of single supplied field
if self.first is self.second:
c2 = c1
else:
c2 = self.second.paint(mode='complex', Nmesh=self.attrs['Nmesh'])
c2 = self.second.compute(mode='complex', Nmesh=self.attrs['Nmesh'])

# calculate the 3d power spectrum, slab-by-slab to save memory
p3d = c1
Expand Down Expand Up @@ -476,7 +476,7 @@ def run(self):
- modes :
the number of Fourier modes averaged together in each bin
"""
c1 = self.first.paint(Nmesh=self.attrs['Nmesh'], mode='complex')
c1 = self.first.compute(Nmesh=self.attrs['Nmesh'], mode='complex')
r1 = c1.preview(self.attrs['Nmesh'], axes=self.attrs['axes'])
# average along projected axes;
# part of product is the rfftn vs r2c (for axes)
Expand All @@ -487,7 +487,7 @@ def run(self):
if self.first is self.second:
c2 = c1
else:
c2 = self.second.paint(Nmesh=self.attrs['Nmesh'], mode='complex')
c2 = self.second.compute(Nmesh=self.attrs['Nmesh'], mode='complex')
r2 = c2.preview(self.attrs['Nmesh'], axes=self.attrs['axes'])
c2 = numpy.fft.rfftn(r2) / self.attrs['Nmesh'].prod() # average along projected axes

Expand Down
4 changes: 2 additions & 2 deletions nbodykit/algorithms/fftrecon.py
Original file line number Diff line number Diff line change
Expand Up @@ -132,7 +132,7 @@ def to_real_field(self):
def run(self):

s_d, s_r = self._compute_s()
return self._paint(s_d, s_r)
return self._helper_paint(s_d, s_r)

def work_with(self, cat, s):
pm = self.pm
Expand Down Expand Up @@ -167,7 +167,7 @@ def _summary_field(self, field, name):
self.logger.info("painted %s, mean=%g" % (name, cmean))


def _paint(self, s_d, s_r):
def _helper_paint(self, s_d, s_r):
""" Convert the displacements of data and random to a single reconstruction mesh object. """

def LGS(delta_s_r):
Expand Down
Loading