Commit e7d7c2c2 authored by PIOLLE's avatar PIOLLE
Browse files

fixes for cerbere v2

parent f5563d99
...@@ -35,7 +35,7 @@ except ImportError: ...@@ -35,7 +35,7 @@ except ImportError:
pyximport.install(setup_args={'include_dirs': numpy.get_include()}) pyximport.install(setup_args={'include_dirs': numpy.get_include()})
# from cerbereutils.resampling import fastregrid # from cerbereutils.resampling import fastregrid
import cerinterp.fastregrid as fastregrid import cerinterp.fastregrid as fastregridR
logging.warning("'pyx' is available but built inplace.") logging.warning("'pyx' is available but built inplace.")
...@@ -407,6 +407,7 @@ def closest_neighbour(source, ...@@ -407,6 +407,7 @@ def closest_neighbour(source,
new=False, new=False,
add_reference=False, add_reference=False,
add_geolocation=True, add_geolocation=True,
add_distance=False,
prefix='resampled_' prefix='resampled_'
): ):
"""Return a new feature corresponding to the resampling of the `source` """Return a new feature corresponding to the resampling of the `source`
...@@ -551,6 +552,11 @@ def closest_neighbour(source, ...@@ -551,6 +552,11 @@ def closest_neighbour(source,
# two-dimensional array # two-dimensional array
else: else:
if numpy.issubdtype(values.dtype, numpy.datetime64):
kd_fillvalue = 0
else:
kd_fillvalue = None
res = kd_tree.get_sample_from_neighbour_info( res = kd_tree.get_sample_from_neighbour_info(
'nn', 'nn',
target_def.shape, target_def.shape,
...@@ -558,15 +564,21 @@ def closest_neighbour(source, ...@@ -558,15 +564,21 @@ def closest_neighbour(source,
valid_input_index, valid_input_index,
valid_output_index, valid_output_index,
index_array, index_array,
fill_value=None fill_value=kd_fillvalue)
)
res = numpy.ma.fix_invalid(res)
if numpy.issubdtype(values.dtype, numpy.datetime64):
deft = numpy.array([0]).astype(res.dtype)
res = numpy.ma.masked_equal(res, deft, copy=False)
target_shape = target.geodimsizes target_shape = target.geodimsizes
if isinstance(target, Grid): if isinstance(target, Grid):
res = res.reshape(target_shape) res = res.reshape(target_shape)
new_dims = target.geodimsizes new_dims = target.geodimsizes
if fieldname in ['time', 'lat', 'lon']: if fieldname in ['time', 'lat', 'lon']:
if not res.any(): if res.count() == 0:
# return None if no valid source data can be resampled on # return None if no valid source data can be resampled on
# the target # the target
if not (values.min() == 0. and values.max() == 0): if not (values.min() == 0. and values.max() == 0):
...@@ -585,16 +597,16 @@ def closest_neighbour(source, ...@@ -585,16 +597,16 @@ def closest_neighbour(source,
field = Field( field = Field(
res, res,
name=newfieldname, name=newfieldname,
dims=new_dims, dims=target.geodimnames,
units=sourcefield.units, units=sourcefield.units,
attrs=copy.copy(sourcefield.attributes), attrs=copy.copy(sourcefield.attrs),
fillvalue=sourcefield.fillvalue) fillvalue=sourcefield.fill_value)
outfeature.add_field(field) outfeature.add_field(field)
if add_reference: if add_reference:
# Add the indices (in the source feature data) of the # Add the indices (in the source feature data) of the
# selected pixels # selected pixels
target_dims = target.geodimsizes target_dims = target.geodimnames
if isinstance(source, Swath) or isinstance(source, Image) or isinstance( if isinstance(source, Swath) or isinstance(source, Image) or isinstance(
source, Trajectory): source, Trajectory):
xdimname = 'cell' xdimname = 'cell'
...@@ -634,7 +646,7 @@ def closest_neighbour(source, ...@@ -634,7 +646,7 @@ def closest_neighbour(source,
fieldy = Field( fieldy = Field(
reference_y, reference_y,
name='resampled_%{}'.format(ydimname), name='resampled_{}'.format(ydimname),
description='{} index of the selected pixel'.format(ydimname), description='{} index of the selected pixel'.format(ydimname),
dims=target_dims, dims=target_dims,
dtype=numpy.dtype(numpy.int32), dtype=numpy.dtype(numpy.int32),
...@@ -643,12 +655,13 @@ def closest_neighbour(source, ...@@ -643,12 +655,13 @@ def closest_neighbour(source,
outfeature.add_field(fieldx) outfeature.add_field(fieldx)
outfeature.add_field(fieldy) outfeature.add_field(fieldy)
if add_distance:
# Add the distance # Add the distance
fielddist = Field( fielddist = Field(
numpy.ma.array(distance_array.reshape(target_shape), mask=res.mask), numpy.ma.array(distance_array.reshape(target_shape), mask=res.mask),
name='distance_to_pixel', name='distance_to_pixel',
description='distance of the selected pixel to the pixel center', description='distance of the selected pixel to the pixel center',
dims=target_dims, dims=target.geodimnames,
units='m', units='m',
dtype=numpy.dtype(numpy.float32) dtype=numpy.dtype(numpy.float32)
) )
......
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