Commit 413b6f84 authored by Jait Dixit's avatar Jait Dixit
Browse files

WIP: Add smoothing functionality to HPSpace and GLSpace

parent 2a8190aa
...@@ -68,7 +68,8 @@ class SmoothOperator(EndomorphicOperator): ...@@ -68,7 +68,8 @@ class SmoothOperator(EndomorphicOperator):
kernel = space_obj.distance_array( kernel = space_obj.distance_array(
x.val.get_axes_local_distribution_strategy(axes=axes)) x.val.get_axes_local_distribution_strategy(axes=axes))
kernel = kernel.apply_scalar_function( kernel = kernel.apply_scalar_function(
space_obj.get_codomain_smoothing_function(self.sigma)) space_obj.codomain_smoothing_function(self.sigma,
transform.target))
# transform # transform
smooth_out = transform(smooth_out, spaces=spaces[0]) smooth_out = transform(smooth_out, spaces=spaces[0])
......
...@@ -3,7 +3,7 @@ from __future__ import division ...@@ -3,7 +3,7 @@ from __future__ import division
import itertools import itertools
import numpy as np import numpy as np
from d2o import STRATEGIES as DISTRIBUTION_STRATEGIES from d2o import distributed_data_object, STRATEGIES as DISTRIBUTION_STRATEGIES
from nifty.spaces.space import Space from nifty.spaces.space import Space
from nifty.config import about, nifty_configuration as gc,\ from nifty.config import about, nifty_configuration as gc,\
...@@ -135,9 +135,8 @@ class GLSpace(Space): ...@@ -135,9 +135,8 @@ class GLSpace(Space):
nlat = self.nlat nlat = self.nlat
weight = np.array(list(itertools.chain.from_iterable( weight = np.array(list(itertools.chain.from_iterable(
itertools.repeat(x ** power, nlon) itertools.repeat(x ** power, nlon)
for x in gl.vol(nlat)) for x in gl.vol(nlat))))
))
if axes is not None: if axes is not None:
# reshape the weight array to match the input shape # reshape the weight array to match the input shape
...@@ -154,6 +153,33 @@ class GLSpace(Space): ...@@ -154,6 +153,33 @@ class GLSpace(Space):
return result_x return result_x
def distance_array(self, distribution_strategy):
shape = self.shape[0]
dists = distributed_data_object(
global_shape=self.shape,
dtype=np.float128,
distribution_strategy=distribution_strategy
)
center_lat = np.random.randint(self.nlat)
center_lon = np.random.randint(self.nlon)
for i in range(self.nlat):
for j in range(self.nlon):
dists[i + j] = np.arccos(np.sin(i) * np.sin(center_lat) +
np.cos(i) * np.cos(center_lat) *
np.cos(j - center_lon))
return dists
def codomain_smoothing_function(self, sigma, target):
if sigma is None:
sigma = np.sqrt(2) * np.pi / (target.lmax + 1)
return lambda x: np.exp(-0.5 * x * (x + 1) * sigma**2)
# ---Added properties and methods--- # ---Added properties and methods---
@property @property
...@@ -184,11 +210,3 @@ class GLSpace(Space): ...@@ -184,11 +210,3 @@ class GLSpace(Space):
"WARNING: nlon was set to an unrecommended value: " "WARNING: nlon was set to an unrecommended value: "
"nlon <> 2*nlat-1.") "nlon <> 2*nlat-1.")
return nlon return nlon
...@@ -35,6 +35,7 @@ from __future__ import division ...@@ -35,6 +35,7 @@ from __future__ import division
import numpy as np import numpy as np
from d2o import distributed_data_object, STRATEGIES as DISTRIBUTION_STRATEGIES
from nifty.config import about from nifty.config import about
from nifty.spaces.space import Space from nifty.spaces.space import Space
from nifty.config import nifty_configuration as gc, \ from nifty.config import nifty_configuration as gc, \
...@@ -125,6 +126,40 @@ class HPSpace(Space): ...@@ -125,6 +126,40 @@ class HPSpace(Space):
self._nside = self._parse_nside(nside) self._nside = self._parse_nside(nside)
def distance_array(self, distribution_strategy):
"""
Calculates distance from center to all the points on the sphere
Parameters
----------
distribution_strategy: Result d2o's distribution strategy
Returns
-------
dists: distributed_data_object
"""
# HPSpace is always 1-dimensional
shape = self.shape[0]
dists = distributed_data_object(
global_shape=shape,
dtype=np.float128,
distribution_strategy=distribution_strategy
)
center_vec = hp.pix2vec(self.nside, np.random.randint(shape))
for i in range(shape):
dists[i] = np.arccos(np.dot(hp.pix2vec(self.nside, i), center_vec))
return dists
def codomain_smoothing_function(self, sigma, target):
if sigma is None:
sigma = np.sqrt(2) * np.pi / (target.lmax + 1)
return lambda x: np.exp(-0.5 * x * (x + 1) * sigma**2)
# ---Mandatory properties and methods--- # ---Mandatory properties and methods---
@property @property
......
...@@ -317,7 +317,7 @@ class RGSpace(Space): ...@@ -317,7 +317,7 @@ class RGSpace(Space):
temp[:] = zerocenter temp[:] = zerocenter
return tuple(temp) return tuple(temp)
def get_codomain_smoothing_function(self, sigma): def codomain_smoothing_function(self, sigma, target):
if sigma is None: if sigma is None:
sigma = np.sqrt(2) * np.max(self.distances) sigma = np.sqrt(2) * np.max(self.distances)
......
...@@ -273,6 +273,11 @@ class Space(object): ...@@ -273,6 +273,11 @@ class Space(object):
raise NotImplementedError(about._errors.cstring( raise NotImplementedError(about._errors.cstring(
"ERROR: There is no generic distance_array for Space base class.")) "ERROR: There is no generic distance_array for Space base class."))
def codomain_smoothing_function(self, sigma, target):
raise NotImplementedError(about._errors.cstring(
"ERROR: There is no generic smoothing fuction for Space base class."
))
def hermitian_decomposition(self, x, axes=None): def hermitian_decomposition(self, x, axes=None):
raise NotImplementedError raise NotImplementedError
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment