Skip to content

Commit

Permalink
error handling :monkax:
Browse files Browse the repository at this point in the history
  • Loading branch information
alexeybelkov committed Dec 24, 2023
1 parent a47c6d3 commit a08a23f
Show file tree
Hide file tree
Showing 5 changed files with 47 additions and 14 deletions.
5 changes: 3 additions & 2 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -59,6 +59,7 @@ Works faster only for `ndim<=4, dtype=float32 or float64 (and bool-int16-32-64 i
```python
from imops import interp1d # same as `scipy.interpolate.interp1d`
```
Works faster only for `ndim<=3, dtype=float32 or float64, order=1`

### Fast 2d linear interpolation
From https://github.com/alexeybelkov/parinterp/tree/rework
Expand All @@ -72,12 +73,12 @@ x_points = points[: n // 2]
values = np.random.uniform(low=0.0, high=1.0, size=(len(x_points),))
interp_points = points[n // 2:]
n_jobs = -1 # will be equal to num of CPU cores
interpolator = Linear2DInterpolator(x_points, values, n_jobs)
# You can optionally pass your own triangulation as an np.array of shape [num_triangles, 3], element at (i, j) position is an index of a point from x_points
interpolator = Linear2DInterpolator(x_points, values, n_jobs=n_jobs, triangles=None)
# Also you can pass values to __call__ and rewrite the ones that were passed to __init__
interp_values = interpolator(interp_points, values + 1.0, fill_value=0.0)
```

Works faster only for `ndim<=3, dtype=float32 or float64, order=1`
### Fast binary morphology

```python
Expand Down
2 changes: 1 addition & 1 deletion imops/cpp/src/interp2d/interpolator.h
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,7 @@ class Interpolator {
throw std::invalid_argument("Length mismatch between int_points and their neighbors");
}

size_t n = int_points.shape()[0];
size_t n = int_points.shape(0);
std::vector<double> int_values(n);

omp_set_dynamic(0); // Explicitly disable dynamic teams
Expand Down
11 changes: 9 additions & 2 deletions imops/cpp/src/interp2d/triangulator.h
Original file line number Diff line number Diff line change
Expand Up @@ -24,10 +24,17 @@ class Triangulator {

Triangulator(const pyarr_size_t& pypoints, int n_jobs,
std::optional<pyarr_size_t> pytriangles) {
if (n_jobs < 0 and n_jobs != -1) {
if (n_jobs <= 0 and n_jobs != -1) {
throw std::invalid_argument(
"Invalid number of workers, have to be -1 or positive integer");
"Invalid number of workers, has to be -1 or positive integer");
}
if (pytriangles.has_value()) {
if (pytriangles->shape(1) != 3 or pytriangles->shape(0) == 0 or
pytriangles->size() / 3 != pytriangles->shape(0)) {
throw std::invalid_argument("Passed triangles argument has an incorrect shape");
}
}

n_jobs_ = n_jobs == -1 ? omp_get_num_procs() : n_jobs;
size_t n = pypoints.shape(0);

Expand Down
31 changes: 28 additions & 3 deletions imops/interp2d.py
Original file line number Diff line number Diff line change
Expand Up @@ -23,15 +23,40 @@ class Linear2DInterpolator(Linear2DInterpolatorCpp):
# """

def __init__(self, points: np.ndarray, values: np.ndarray, n_jobs: int = 1, triangles: Optional[np.ndarray] = None, **kwargs):

if triangles is not None:
if isinstance(triangles, np.ndarray):
if triangles.ndim != 2 or triangles.shape[1] != 3 or triangles.shape[0] * 3 != triangles.size:
raise ValueError('Passed \"triangles\" argument has an incorrect shape')
else:
raise ValueError(f'Wrong type of \"triangles\" argument, expected Optional[np.ndarray]. Got {type(triangles)}')

if n_jobs <= 0 and n_jobs != -1:
raise ValueError('Invalid number of workers, has to be -1 or positive integer')

super().__init__(points, n_jobs, triangles)
self.kdtree = KDTree(data=points, **kwargs)
self.values = values
self.n_jobs = n_jobs

def __call__(
self, points: np.ndarray, values: Optional[np.ndarray] = None, fill_value: float = 0.0
) -> np.ndarray[np.float32]:
def __call__(self, points: np.ndarray, values: Optional[np.ndarray] = None, fill_value: float = 0.0) -> np.ndarray[np.float32]:

if values is not None:
if isinstance(values, np.ndarray):
if values.ndim > 1:
raise ValueError(f'Wrong shape of \"values\" argument, expected ndim=1. Got shape {values.shape}')
else:
raise ValueError(f'Wrong type of \"values\" argument, expected Optional[np.ndarray]. Got {type(values)}')

self.values = self.values if values is None else values

if self.values is None:
raise ValueError('\"values\" argument was never passed neither in __init__ or __call__ methods')

_, neighbors = self.kdtree.query(points, 1, workers=self.n_jobs)

if not isinstance(points, np.ndarray):
raise ValueError(f'Wrong type of \"points\" argument, expected np.ndarray. Got {type(points)}')

int_values = super().__call__(points, self.values, neighbors, fill_value)
return int_values
12 changes: 6 additions & 6 deletions tests/test_interp2d.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,33 +6,33 @@
from deli import load_numpy

np.random.seed(1337)

n_jobs = 1

def test_test_data():
for i in range(2):
x = load_numpy(f'test_data/arr_{i}.npy')
int_points = np.transpose((np.isnan(x)).nonzero())
x_points: np.ndarray = np.transpose((~np.isnan(x)).nonzero())
x_values = x[~np.isnan(x)]
imops_values = Linear2DInterpolator(x_points, x_values, n_jobs=1)(int_points, fill_value=0.0)
imops_values = Linear2DInterpolator(x_points, x_values, n_jobs=n_jobs)(int_points, fill_value=0.0)
triangles = scipy.spatial.Delaunay(x_points).simplices
delaunay_values = Linear2DInterpolator(x_points, x_values, n_jobs=1, triangles=triangles)(int_points, fill_value=0.0)
delaunay_values = Linear2DInterpolator(x_points, x_values, n_jobs=n_jobs, triangles=triangles)(int_points, fill_value=0.0)
scipy_values = griddata(x_points, x_values, int_points, method='linear', fill_value=0.0)

delta_ds = np.abs(delaunay_values - scipy_values)
# delta_di = np.abs(delaunay_values - imops_values)
delta_si = np.abs(scipy_values - imops_values)

assert delta_ds.max() <= 1e-8 and delta_si.max() <= 5, 'Failed with big cases, arr_0 or arr_1'
assert delta_ds.max() <= 1e-10 and delta_si.max() <= 5, 'Failed with big cases, arr_0 or arr_1'

for i in range(2, 7):
x = load_numpy(f'test_data/arr_{i}.npy')
x_values = load_numpy(f'test_data/val_{i}.npy')
x_points = np.transpose(x.nonzero())
int_points = np.transpose((~x).nonzero())
imops_values = Linear2DInterpolator(x_points, x_values, n_jobs=1)(int_points, fill_value=0.0)
imops_values = Linear2DInterpolator(x_points, x_values, n_jobs=n_jobs)(int_points, fill_value=0.0)
triangles = scipy.spatial.Delaunay(x_points).simplices
delaunay_values = Linear2DInterpolator(x_points, x_values, n_jobs=1, triangles=triangles)(int_points, fill_value=0.0)
delaunay_values = Linear2DInterpolator(x_points, x_values, n_jobs=n_jobs, triangles=triangles)(int_points, fill_value=0.0)
scipy_values = griddata(x_points, x_values, int_points, method='linear', fill_value=0.0)

delta_ds = np.abs(delaunay_values - scipy_values)
Expand Down

0 comments on commit a08a23f

Please sign in to comment.