gltransform.py 2.22 KB
Newer Older
1
2
3
4
5
6
7
8
import numpy as np
from transform import Transform
from d2o import distributed_data_object
from nifty.config import dependency_injector as gdi
import nifty.nifty_utilities as utilities

gl = gdi.get('libsharp_wrapper_gl')

9

10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
class GLTransform(Transform):
    """
        GLTransform wrapper for libsharp's transform functions
    """

    def __init__(self, domain, codomain):
        self.domain = domain
        self.codomain = codomain

        if 'libsharp_wrapper_gl' not in gdi:
            raise ImportError("The module libsharp_wrapper_gl " +
                              "is needed but not available")

    def transform(self, val, axes, **kwargs):
        if self.domain.discrete:
25
            val = self.domain.calc_weight(val, power=-0.5)
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67

        # shorthands for transform parameters
        nlat = self.domain.paradict['nlat']
        nlon = self.domain.paradict['nlon']
        lmax = self.codomain.paradict['lmax']
        mmax = self.codomain.paradict['mmax']

        if isinstance(val, distributed_data_object):
            temp_val = val.get_full_data()
        else:
            temp_val = val

        return_val = None

        for slice_list in utilities.get_slice_list(temp_val.shape, axes):
            if slice_list == [slice(None, None)]:
                inp = temp_val
            else:
                if return_val is None:
                    return_val = np.empty_like(temp_val)
                inp = temp_val[slice_list]

            if self.domain.dtype == np.dtype('float32'):
                inp = gl.map2alm_f(inp,
                                   nlat=nlat, nlon=nlon,
                                   lmax=lmax, mmax=mmax)
            else:
                inp = gl.map2alm(inp,
                                 nlat=nlat, nlon=nlon,
                                 lmax=lmax, mmax=mmax)

            if slice_list == [slice(None, None)]:
                return_val = inp
            else:
                return_val[slice_list] = inp

        if isinstance(val, distributed_data_object):
            new_val = val.copy_empty(dtype=self.codomain.dtype)
            new_val.set_full_data(return_val, copy=False)
        else:
            return_val = return_val.astype(self.codomain.dtype, copy=False)

68
        return return_val