Skip to content

Commit

Permalink
Merge pull request #44 from EMeraldGeo/optimization
Browse files Browse the repository at this point in the history
Optimizations that do NOT change functionality
  • Loading branch information
mmaelicke authored Sep 18, 2020
2 parents 5ef4696 + 8031dc3 commit e09a650
Show file tree
Hide file tree
Showing 2 changed files with 53 additions and 60 deletions.
78 changes: 43 additions & 35 deletions skgstat/DirectionalVariogram.py
Original file line number Diff line number Diff line change
Expand Up @@ -194,6 +194,11 @@ def __init__(self,
Set the Verbosity of the class. Not Implemented yet.
"""

# FIXME: Call __init__ of baseclass?

self._direction_mask_cache = None

# Set coordinates
self._X = np.asarray(coordinates)

Expand All @@ -205,7 +210,8 @@ def __init__(self,

# set values
self._values = None
self.set_values(values=values)
# calc_diff = False here, because it will be calculated by fit() later
self.set_values(values=values, calc_diff=False)

# distance matrix
self._dist = None
Expand Down Expand Up @@ -268,7 +274,7 @@ def __init__(self,
self._cache_experimental = False

# do the preprocessing and fitting upon initialization
self.preprocessing(force=True)
# Note that fit() calls preprocessing
self.fit(force=True)

@property
Expand Down Expand Up @@ -452,20 +458,21 @@ def local_reference_system(self, poi):
the original coordinates.
"""
# define a point-wise transform function
def _transform(p1, p2, a):
p = p1 - p2
x = p[0] * np.cos(a) - p[1] * np.sin(a)
y = p[0] * np.sin(a) + p[1] * np.cos(a)
return np.array([x, y])



# get the azimuth in radians
gamma = np.radians(self.azimuth)

# transform
_X = np.fromiter(chain.from_iterable(
map(lambda p: _transform(p, poi, gamma), self._X)
), dtype=self._X.dtype).reshape(self._X.shape)
cosgamma = np.cos(gamma)
singamma = np.sin(gamma)

_X = np.zeros(dtype=self._X.dtype, shape=self._X.shape)

for i in range(0, len(self._X)):
p = self._X[i,:] - poi
_X[i,0] = p[0] * cosgamma - p[1] * singamma
_X[i,1] = p[0] * singamma + p[1] * cosgamma

# return
return _X
Expand All @@ -488,17 +495,14 @@ def _calc_groups(self, force=False):
self._groups[np.where(~self._direction_mask())] = -1

# @jit
def _direction_mask(self):
def _direction_mask(self, force=False):
"""Directional Mask
Array aligned to self.distance masking all point pairs which shall be
ignored for binning and grouping. The one dimensional array contains
all row-wise point pair combinations from the upper or lower triangle
of the distance matrix in case either of both is directional.
TODO: This array is not cached. it is used twice, for binning and
grouping.
Returns
-------
mask : numpy.array
Expand All @@ -507,30 +511,34 @@ def _direction_mask(self):
not.
"""
# build the full coordinate matrix
n = len(self._X)
_mask = np.zeros((n, n), dtype=bool)

# build the masking
for i in range(n):
loc = self.local_reference_system(poi=self._X[i])
if force or self._direction_mask_cache is None:
# build the full coordinate matrix
n = len(self._X)
_mask = np.zeros((n, n), dtype=bool)

# apply the search radius
sr = self._directional_model(local_ref=loc)
# build the masking
for i in range(n):
loc = self.local_reference_system(poi=self._X[i])

_m = np.fromiter(
(Point(p).within(sr) or Point(p).touches(sr) for p in loc),
dtype=bool)
_mask[:, i] = _m
# apply the search radius
sr = self._directional_model(local_ref=loc)

# combine lower and upper triangle
def _indexer():
for i in range(n):
for j in range(n):
if i < j:
yield _mask[i, j] or _mask[j, i]
_m = np.fromiter(
(Point(p).within(sr) or Point(p).touches(sr) for p in loc),
dtype=bool)
_mask[:, i] = _m

# combine lower and upper triangle
def _indexer():
for i in range(n):
for j in range(n):
if i < j:
yield _mask[i, j] or _mask[j, i]

self._direction_mask_cache = np.fromiter(_indexer(), dtype=bool)

return np.fromiter(_indexer(), dtype=bool)
return self._direction_mask_cache

def search_area(self, poi=0, ax=None):
"""Plot Search Area
Expand Down
35 changes: 10 additions & 25 deletions skgstat/Variogram.py
Original file line number Diff line number Diff line change
Expand Up @@ -162,7 +162,8 @@ def __init__(self,

# set values
self._values = None
self.set_values(values=values)
# calc_diff = False here, because it will be calculated by fit() later
self.set_values(values=values, calc_diff=False)

# distance matrix
self._dist = None
Expand Down Expand Up @@ -214,7 +215,7 @@ def __init__(self,
self._cache_experimental = False

# do the preprocessing and fitting upon initialization
self.preprocessing(force=True)
# Note that fit() calls preprocessing
self.fit(force=True)

@property
Expand Down Expand Up @@ -276,7 +277,7 @@ def value_matrix(self):
"""
return squareform(self._diff)

def set_values(self, values):
def set_values(self, values, calc_diff=True):
"""Set new values
Will set the passed array as new value array. This array has to be of
Expand Down Expand Up @@ -328,7 +329,8 @@ class does only accept one dimensional arrays.
self._values = np.asarray(values)

# recalculate the pairwise differences
self._calc_diff(force=True)
if calc_diff:
self._calc_diff(force=True)

@property
def bin_func(self):
Expand Down Expand Up @@ -1039,27 +1041,11 @@ def _calc_diff(self, force=False):
self._diff = np.zeros(int((l**2 - l) / 2))

# calculate the pairwise differences
for t, k in zip(self.__vdiff_indexer(), range(len(self._diff))):
self._diff[k] = np.abs(v[t[0]] - v[t[1]])

#@jit
def __vdiff_indexer(self):
"""Pairwise indexer
Returns an iterator over the values or coordinates in squareform
coordinates. The iterable will be of type tuple.
Returns
-------
iterable
"""
l = len(self.values)

k = 0
for i in range(l):
for j in range(l):
if i < j:
yield i, j
for j in range(i+1, l):
self._diff[k] = np.abs(v[i] - v[j])
k += 1

def _calc_groups(self, force=False):
"""Calculate the lag class mask array
Expand Down Expand Up @@ -1119,7 +1105,6 @@ def experimental(self):
return self._experimental

@property
# @jit
def _experimental(self):
"""
Expand Down

0 comments on commit e09a650

Please sign in to comment.