Commit 1e6af087 authored by Ultima's avatar Ultima
Browse files

Cleaning up. Making more files PEP8 compliant.

Abstracted the power_indices class from rg_space -> Implemented a lm_power_indices class.
In the middle of space-interface cleaning.
parent 6c50de1e
......@@ -35,7 +35,6 @@ from keepers import about,\
from nifty_cmaps import ncmap
from nifty_core import space,\
point_space,\
nested_space,\
field
from nifty_mpi_data import distributed_data_object, d2o_librarian
......
......@@ -16,7 +16,11 @@ MAX = np.max
SUM = np.sum
class _COMM_WORLD():
class Comm(object):
pass
class Intracomm(Comm):
def __init__(self):
self.rank = 0
self.size = 1
......@@ -111,4 +115,4 @@ COMPLEX = _datatype("MPI_COMPLEX")
DOUBLE_COMPLEX = _datatype("MPI_DOUBLE_COMPLEX")
COMM_WORLD = _COMM_WORLD()
COMM_WORLD = Intracomm()
# -*- coding: utf-8 -*-
import ConfigParser
import pprint
class variable(object):
def __init__(self, name, default_value_list, checker, genus='str'):
......@@ -47,7 +49,7 @@ class variable(object):
return self.name
def __repr__(self):
return "<nifty variable '" + str(self.name) + "': " + \
return "<" + str(self.name) + "': " + \
str(self.get_value()) + ">"
......@@ -161,7 +163,9 @@ class configuration(object):
item.set_value(temp_value)
def __repr__(self):
return "<nifty configuration> \n" + self.variable_dict.__repr__()
return_string = ("<nifty configuration> \n" +
pprint.pformat(self.variable_dict))
return return_string
......
......@@ -6,6 +6,9 @@ from nifty_dependency_injector import dependency_injector
from nifty_configuration import variable,\
configuration
global_dependency_injector = dependency_injector(
['h5py',
('mpi4py.MPI', 'MPI'),
......@@ -23,12 +26,12 @@ variable_fft_module = variable('fft_module',
variable_lm2gl = variable('lm2gl',
[True, False],
lambda z: z is True or z is False,
'boolean')
genus = 'boolean')
variable_verbosity = variable('verbosity',
[1],
lambda z: z == abs(int(z)),
'int')
genus = 'int')
variable_mpi_module = variable('mpi_module',
['MPI', 'MPI_dummy'],
......@@ -49,3 +52,12 @@ global_configuration = configuration(
variable_default_distribution_strategy
],
path=os.path.expanduser('~') + "/.nifty/global_config")
variable_default_comm = variable(
'default_comm',
['COMM_WORLD'],
lambda z: hasattr(global_dependency_injector[
global_configuration['mpi_module']], z))
global_configuration.register(variable_default_comm)
## NIFTY (Numerical Information Field Theory) has been developed at the
## Max-Planck-Institute for Astrophysics.
##
## Copyright (C) 2015 Max-Planck-Society
##
## Author: Marco Selig
## Project homepage: <http://www.mpa-garching.mpg.de/ift/nifty/>
##
## This program is free software: you can redistribute it and/or modify
## it under the terms of the GNU General Public License as published by
## the Free Software Foundation, either version 3 of the License, or
## (at your option) any later version.
##
## This program is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
## See the GNU General Public License for more details.
##
## You should have received a copy of the GNU General Public License
## along with this program. If not, see <http://www.gnu.org/licenses/>.
# NIFTY (Numerical Information Field Theory) has been developed at the
# Max-Planck-Institute for Astrophysics.
#
# Copyright (C) 2015 Max-Planck-Society
#
# Author: Marco Selig
# Project homepage: <http://www.mpa-garching.mpg.de/ift/nifty/>
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
# See the GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program. If not, see <http://www.gnu.org/licenses/>.
"""
.. __ ____ __
......@@ -32,43 +32,36 @@
"""
from __future__ import division
#from nifty import *
import os
import numpy as np
from numpy import pi
import pylab as pl
from matplotlib.colors import LogNorm as ln
from matplotlib.ticker import LogFormatter as lf
from nifty.keepers import about
from nifty.nifty_core import pi, \
space, \
point_space, \
field
from nifty.nifty_core import space,\
point_space,\
field
from keepers import about,\
global_configuration as gc,\
global_dependency_injector as gdi
from nifty.nifty_paradict import lm_space_paradict,\
gl_space_paradict,\
hp_space_paradict
gl_space_paradict,\
hp_space_paradict
from nifty.nifty_power_indices import lm_power_indices
from nifty.nifty_random import random
#import libsharp_wrapper_gl as gl
try:
import libsharp_wrapper_gl as gl
except(ImportError):
about.infos.cprint("INFO: libsharp_wrapper_gl not available.")
_gl_available = False
gl = gdi['libsharp_wrapper_gl']
hp = gdi['healpy']
if gl is None and gc['lm2gl']:
about.warnings.cprint("WARNING: global setting 'about.lm2gl' corrected.")
about.lm2gl.off()
else:
_gl_available = True
#import healpy as hp
try:
import healpy as hp
except(ImportError):
about.infos.cprint("INFO: healpy not available.")
_hp_available = False
else:
_hp_available = True
gc['lm2gl'] = False
LM_DISTRIBUTION_STRATEGIES = []
##-----------------------------------------------------------------------------
class lm_space(point_space):
"""
......@@ -90,7 +83,7 @@ class lm_space(point_space):
mmax : int, *optional*
Maximum :math:`m`-value up to which the spherical harmonics
coefficients are to be used (default: `lmax`).
datatype : numpy.dtype, *optional*
dtype : numpy.dtype, *optional*
Data type of the field values (default: numpy.complex128).
See Also
......@@ -120,7 +113,7 @@ class lm_space(point_space):
para : numpy.ndarray
One-dimensional array containing the two numbers `lmax` and
`mmax`.
datatype : numpy.dtype
dtype : numpy.dtype
Data type of the field values.
discrete : bool
Parameter captioning the fact that an :py:class:`lm_space` is
......@@ -128,7 +121,9 @@ class lm_space(point_space):
vol : numpy.ndarray
Pixel volume of the :py:class:`lm_space`, which is always 1.
"""
def __init__(self, lmax, mmax=None, datatype=None, datamodel = 'np'):
def __init__(self, lmax, mmax=None, dtype=np.dtype('complex128'),
datamodel='np'):
"""
Sets the attributes for an lm_space class instance.
......@@ -140,7 +135,7 @@ class lm_space(point_space):
mmax : int, *optional*
Maximum :math:`m`-value up to which the spherical harmonics
coefficients are to be used (default: `lmax`).
datatype : numpy.dtype, *optional*
dtype : numpy.dtype, *optional*
Data type of the field values (default: numpy.complex128).
Returns
......@@ -157,21 +152,21 @@ class lm_space(point_space):
"""
## check imports
if(not _gl_available)and(not _hp_available):
raise ImportError(about._errors.cstring("ERROR: neither libsharp_wrapper_gl nor healpy available."))
# check imports
if 'libsharp_wrapper_gl' not in gdi and 'healpy' not in gdi:
raise ImportError(about._errors.cstring(
"ERROR: neither libsharp_wrapper_gl nor healpy available."))
self.paradict = lm_space_paradict(lmax=lmax, mmax=mmax)
## check data type
if(datatype is None):
datatype = np.complex128
elif(datatype not in [np.complex64,np.complex128]):
# check data type
dtype = np.dtype(dtype)
if dtype not in [np.dtype('complex64'), np.dtype('complex128')]:
about.warnings.cprint("WARNING: data type set to default.")
datatype = np.complex128
self.datatype = datatype
dtype = np.dtype('complex128')
self.dtype = dtype
## set datamodel
# set datamodel
if datamodel not in ['np']:
about.warnings.cprint("WARNING: datamodel set to default.")
self.datamodel = 'np'
......@@ -180,7 +175,14 @@ class lm_space(point_space):
self.discrete = True
self.harmonic = True
self.vol = np.real(np.array([1],dtype=self.datatype))
self.distances = (np.float(1),)
self.power_indices = lm_power_indices(
lmax=self.paradict['lmax'],
dim=self.get_dim(),
comm=self.comm,
datamodel=self.datamodel,
allowed_distribution_strategies=LM_DISTRIBUTION_STRATEGIES)
@property
def para(self):
......@@ -188,80 +190,104 @@ class lm_space(point_space):
self.paradict['mmax']], dtype=int)
return temp
@para.setter
def para(self, x):
self.paradict['lmax'] = x[0]
self.paradict['mmax'] = x[1]
def copy(self):
return lm_space(lmax = self.paradict['lmax'],
mmax = self.paradict['mmax'],
datatype = self.datatype)
##+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
def lmax(self):
"""
Returns the maximum quantum number :math:`\ell`.
Returns
-------
lmax : int
Maximum quantum number :math:`\ell`.
"""
return self.paradict['lmax']
def mmax(self):
"""
Returns the maximum quantum number :math:`m`.
Returns
-------
mmax : int
Maximum quantum number :math:`m`.
"""
return self.paradict['mmax']
return lm_space(lmax=self.paradict['lmax'],
mmax=self.paradict['mmax'],
dtype=self.dtype)
# def _enforce_shape(self, x):
# """
# Shapes an array of valid field values correctly, according to the
# specifications of the space instance.
#
# Parameters
# ----------
# x : numpy.ndarray
# Array containing the field values to be put into shape.
#
# Returns
# -------
# y : numpy.ndarray
# Correctly shaped array.
# """
# about.warnings.cprint("WARNING: enforce_shape is deprecated!")
#
# x = np.array(x)
# if(np.size(x) != self.get_dim(split=False)):
# raise ValueError(about._errors.cstring("ERROR: dimension mismatch ( " +
# str(np.size(x)) + " <> " + str(self.get_dim(split=False)) + " )."))
# elif(not np.all(np.array(np.shape(x))==self.get_dim(split=True))):
# about.warnings.cprint("WARNING: reshaping forced.")
#
# return x.reshape(self.get_dim(split=True), order='C')
# def lmax(self):
# """
# Returns the maximum quantum number :math:`\ell`.
#
# Returns
# -------
# lmax : int
# Maximum quantum number :math:`\ell`.
# """
# return self.paradict['lmax']
#
# def mmax(self):
# """
# Returns the maximum quantum number :math:`m`.
#
# Returns
# -------
# mmax : int
# Maximum quantum number :math:`m`.
#
# """
# return self.paradict['mmax']
def get_shape(self):
mmax = self.paradict['mmax']
lmax = self.paradict['lmax']
return np.array([(mmax+1)*(lmax+1)-((lmax+1)*lmax)//2], dtype=int)
def get_dim(self, split=False):
"""
Computes the dimension of the space, i.e.\ the number of spherical
harmonics components that are stored.
Parameters
----------
split : bool, *optional*
Whether to return the dimension as an array with one component
or as a scalar (default: False).
Returns
-------
dim : {int, numpy.ndarray}
Number of spherical harmonics components.
Notes
-----
Due to the symmetry assumption, only the components with
non-negative :math:`m` are stored and only these components are
counted here.
"""
## dim = (mmax+1)*(lmax-mmax/2+1)
if(split):
return self.get_shape()
#return np.array([(self.para[0]+1)*(self.para[1]+1)-(self.para[1]+1)*self.para[1]//2],dtype=np.int)
else:
return np.prod(self.get_shape())
#return (self.para[0]+1)*(self.para[1]+1)-(self.para[1]+1)*self.para[1]//2
##+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
def get_dof(self):
return (np.int((mmax + 1) * (lmax + 1) - ((lmax + 1) * lmax) // 2),)
# -> inherited
# def get_dim(self, split=False):
# """
# Computes the dimension of the space, i.e.\ the number of spherical
# harmonics components that are stored.
#
# Parameters
# ----------
# split : bool, *optional*
# Whether to return the dimension as an array with one component
# or as a scalar (default: False).
#
# Returns
# -------
# dim : {int, numpy.ndarray}
# Number of spherical harmonics components.
#
# Notes
# -----
# Due to the symmetry assumption, only the components with
# non-negative :math:`m` are stored and only these components are
# counted here.
# """
# # dim = (mmax+1)*(lmax-mmax/2+1)
# if(split):
# return self.get_shape()
# # return
# # np.array([(self.para[0]+1)*(self.para[1]+1)-(self.para[1]+1)*self.para[1]//2],dtype=np.int)
# else:
# return np.prod(self.get_shape())
# # return
# # (self.para[0]+1)*(self.para[1]+1)-(self.para[1]+1)*self.para[1]//2
def get_dof(self, split=False):
"""
Computes the number of degrees of freedom of the space, taking into
account symmetry constraints and complex-valuedness.
......@@ -276,107 +302,66 @@ class lm_space(point_space):
The number of degrees of freedom is reduced due to the hermitian
symmetry, which is assumed for the spherical harmonics components.
"""
## dof = 2*dim-(lmax+1) = (lmax+1)*(2*mmax+1)*(mmax+1)*mmax
return (self.para[0]+1)*(2*self.para[1]+1)-(self.para[1]+1)*self.para[1]
##+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
# dof = 2*dim-(lmax+1) = (lmax+1)*(2*mmax+1)*(mmax+1)*mmax
mmax = self.paradict['lmax']
lmax = self.paradict['lmax']
dof = np.int((lmax + 1) * (2 * mmax + 1) - (mmax + 1) * mmax)
if split:
return (dof, )
else:
return dof
def enforce_power(self,spec,**kwargs):
def get_meta_volume(self, split=False):
"""
Provides a valid power spectrum array from a given object.
Parameters
----------
spec : {float, list, numpy.ndarray, nifty.field, function}
Fiducial power spectrum from which a valid power spectrum is to
be calculated. Scalars are interpreted as constant power
spectra.
Returns
-------
spec : numpy.ndarray
Valid power spectrum.
"""
if(isinstance(spec,field)):
spec = spec.val.astype(dtype=self.datatype)
elif(callable(spec)):
try:
spec = np.array(spec(np.arange(self.para[0]+1,dtype=self.vol.dtype)),dtype=self.datatype) ## prevent integer division
except:
raise TypeError(about._errors.cstring("ERROR: invalid power spectra function.")) ## exception in ``spec(kindex)``
elif(np.isscalar(spec)):
spec = np.array([spec],dtype=self.datatype)
else:
spec = np.array(spec,dtype=self.datatype)
Calculates the meta volumes.
## drop imaginary part
spec = np.real(spec)
## check finiteness
if(not np.all(np.isfinite(spec))):
about.warnings.cprint("WARNING: infinite value(s).")
## check positivity (excluding null)
if(np.any(spec<0)):
raise ValueError(about._errors.cstring("ERROR: nonpositive value(s)."))
elif(np.any(spec==0)):
about.warnings.cprint("WARNING: nonpositive value(s).")
size = self.para[0]+1 ## lmax+1
## extend
if(np.size(spec)==1):
spec = spec*np.ones(size,dtype=spec.dtype,order='C')
## size check
elif(np.size(spec)<size):
raise ValueError(about._errors.cstring("ERROR: size mismatch ( "+str(np.size(spec))+" < "+str(size)+" )."))
elif(np.size(spec)>size):
about.warnings.cprint("WARNING: power spectrum cut to size ( == "+str(size)+" ).")
spec = spec[:size]
return spec
##+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
def _getlm(self): ## > compute all (l,m)
index = np.arange(self.get_dim(split=False))
n = 2*self.para[0]+1
m = np.ceil((n-np.sqrt(n**2-8*(index-self.para[0])))/2).astype(np.int)
l = index-self.para[0]*m+m*(m-1)//2
return l,m
def set_power_indices(self,**kwargs):
"""
Sets the (un)indexing objects for spectral indexing internally.
The meta volumes are the volumes associated with each component of
a field, taking into account field components that are not
explicitly included in the array of field values but are determined
by symmetry conditions.
Parameters
----------
None
total : bool, *optional*
Whether to return the total meta volume of the space or the
individual ones of each field component (default: False).
Returns
-------
None
See Also
--------
get_power_indices
mol : {numpy.ndarray, float}
Meta volume of the field components or the complete space.
Notes
-----
The spherical harmonics components with :math:`m=0` have meta
volume 1, the ones with :math:`m>0` have meta volume 2, sinnce they
each determine another component with negative :math:`m`.
"""
## check storage
if(not hasattr(self,"power_indices")):
## power indices
# about.infos.cflush("INFO: setting power indices ...")
kindex = np.arange(self.para[0]+1,dtype=np.int)
rho = 2*kindex+1
if(_hp_available): ## default
pindex = hp.Alm.getlm(self.para[0],i=None)[0] ## l of (l,m)
else:
pindex = self._getlm()[0] ## l of (l,m)
pundex = np.unique(pindex,return_index=True,return_inverse=False)[1]
## storage
self.power_indices = {"kindex":kindex,"pindex":pindex,"pundex":pundex,"rho":rho} ## alphabetical
# about.infos.cprint(" done.")
##+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
if not split:
return np.float(self.get_dof())
else:
mol = self.cast(1, dtype=np.float)
mol[self.paradict['lmax'] + 1:] = 2 # redundant: (l,m) and (l,-m)
return mol
def enforce_values(self,x,extend=True):
# TODO: Extend to binning/log
def enforce_power(self, spec):
size = self.paradict['lmax'] + 1
kindex = np.arange(size, dtype=np.array(self.distances).dtype)
return self._enforce_power_helper(spec=spec,
size=size,
kindex=kindex)
def _getlm(self): # > compute all (l,m)
index = np.arange(self.get_dim())
n = 2 * self.paradict['lmax'] + 1
m = np.ceil(
(n - np.sqrt(n**2 - 8 * (index - self.paradict['lmax']))) / 2
).astype(np.int)
l = index - self.paradict['lmax'] * m + m * (m - 1) // 2
return l, m
def _enforce_values(self, x, extend=True):
"""
Computes valid field values from a given object, taking into
account data types, size, and hermitian symmetry.
......@@ -397,40 +382,42 @@ class lm_space(point_space):
Whether a scalar is extented to a constant array or not
(default: True).
"""
if(isinstance(x,field)):
if(self==x.domain):
if(self.datatype is not x.domain.datatype):
raise TypeError(about._errors.cstring("ERROR: inequal data types ( '"+str(np.result_type(self.datatype))+"' <> '"+str(np.result_type(x.domain.datatype))+"' )."))
if(isinstance(x, field)):
if(self == x.domain):
if(self.dtype is not x.domain.dtype):
raise TypeError(about._errors.cstring("ERROR: inequal data types ( '" + str(
np.result_type(self.dtype)) + "' <> '" + str(np.result_type(x.domain.dtype)) + "' )."))
else:
x = np.copy(x.val)
else:
raise ValueError(about._errors.cstring("ERROR: inequal domains."))
raise ValueError(about._errors.cstring(
"ERROR: inequal domains."))
else:
if(np.size(x)==1):
if(np.size(x) == 1):