Add support in matplotlib.pyplot.imshow for float16

[geometrikal]

Bug report

Bug summary

imshow of array with dtype float16 fails with unsupported dtype

Code for reproduction

import numpy as np
import matplotlib.pyplot as plt
#This works:
a = np.asarray([[0.1,0.2],[0.3,0.4]], dtype=np.float64)
plt.imshow(a), plt.show()
#This fails
a = np.asarray([[0.1,0.2],[0.3,0.4]], dtype=np.float16)
plt.imshow(a), plt.show()

Matplotlib version

• Operating system: Windows 10 1903
• Matplotlib version: 3.1.1
• Matplotlib backend (print(matplotlib.get_backend())): module://backend_interagg
• Python version: 3.73
• Jupyter version (if applicable): n/a
• Other libraries:

installed from default conda

[tacaswell]

In [1]: 
In [2]: import matplotlib

In [3]: matplotlib.use('agg')

In [4]: import matplotlib.pyplot as plt

In [5]: a = np.asarray([[0.1,0.2],[0.3,0.4]], dtype=np.float16)
---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
<ipython-input-5-aa5565166ec5> in <module>()
----> 1 a = np.asarray([[0.1,0.2],[0.3,0.4]], dtype=np.float16)

NameError: name 'np' is not defined

In [6]: import numpy as np

In [7]: a = np.asarray([[0.1,0.2],[0.3,0.4]], dtype=np.float16)

In [8]: im = plt.imshow(a)

In [9]: im.axes.figure.canvas.draw()
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-9-366a74871c20> in <module>()
----> 1 im.axes.figure.canvas.draw()

~/src/p/matplotlib/lib/matplotlib/backends/backend_agg.py in draw(self)
    391              (self.toolbar._wait_cursor_for_draw_cm() if self.toolbar
    392               else nullcontext()):
--> 393             self.figure.draw(self.renderer)
    394             # A GUI class may be need to update a window using this draw, so
    395             # don't forget to call the superclass.

~/src/p/matplotlib/lib/matplotlib/artist.py in draw_wrapper(artist, renderer, *args, **kwargs)
     36                 renderer.start_filter()
     37 
---> 38             return draw(artist, renderer, *args, **kwargs)
     39         finally:
     40             if artist.get_agg_filter() is not None:

~/src/p/matplotlib/lib/matplotlib/figure.py in draw(self, renderer)
   1698             self.patch.draw(renderer)
   1699             mimage._draw_list_compositing_images(
-> 1700                 renderer, self, artists, self.suppressComposite)
   1701 
   1702             renderer.close_group('figure')

~/src/p/matplotlib/lib/matplotlib/image.py in _draw_list_compositing_images(renderer, parent, artists, suppress_composite)
    133     if not_composite or not has_images:
    134         for a in artists:
--> 135             a.draw(renderer)
    136     else:
    137         # Composite any adjacent images together

~/src/p/matplotlib/lib/matplotlib/artist.py in draw_wrapper(artist, renderer, *args, **kwargs)
     36                 renderer.start_filter()
     37 
---> 38             return draw(artist, renderer, *args, **kwargs)
     39         finally:
     40             if artist.get_agg_filter() is not None:

~/src/p/matplotlib/lib/matplotlib/axes/_base.py in draw(self, renderer, inframe)
   2691             renderer.stop_rasterizing()
   2692 
-> 2693         mimage._draw_list_compositing_images(renderer, self, artists)
   2694 
   2695         renderer.close_group('axes')

~/src/p/matplotlib/lib/matplotlib/image.py in _draw_list_compositing_images(renderer, parent, artists, suppress_composite)
    133     if not_composite or not has_images:
    134         for a in artists:
--> 135             a.draw(renderer)
    136     else:
    137         # Composite any adjacent images together

~/src/p/matplotlib/lib/matplotlib/artist.py in draw_wrapper(artist, renderer, *args, **kwargs)
     36                 renderer.start_filter()
     37 
---> 38             return draw(artist, renderer, *args, **kwargs)
     39         finally:
     40             if artist.get_agg_filter() is not None:

~/src/p/matplotlib/lib/matplotlib/image.py in draw(self, renderer, *args, **kwargs)
    623         else:
    624             im, l, b, trans = self.make_image(
--> 625                 renderer, renderer.get_image_magnification())
    626             if im is not None:
    627                 renderer.draw_image(gc, l, b, im)

~/src/p/matplotlib/lib/matplotlib/image.py in make_image(self, renderer, magnification, unsampled)
    880             self._A, bbox, transformed_bbox,
    881             self.get_clip_box() or self.axes.bbox,
--> 882             magnification, unsampled=unsampled)
    883 
    884     def _check_unsampled_image(self, renderer):

~/src/p/matplotlib/lib/matplotlib/image.py in _make_image(self, A, in_bbox, out_bbox, clip_bbox, magnification, unsampled, round_to_pixel_border)
    475                 A_scaled += 0.1
    476                 # resample the input data to the correct resolution and shape
--> 477                 A_resampled = _resample(self, A_scaled, out_shape, t)
    478                 # done with A_scaled now, remove from namespace to be sure!
    479                 del A_scaled

~/src/p/matplotlib/lib/matplotlib/image.py in _resample(image_obj, data, out_shape, transform, resample, alpha)
    200                     alpha,
    201                     image_obj.get_filternorm(),
--> 202                     image_obj.get_filterrad())
    203     return out
    204 

ValueError: Unsupported dtype

The issue is that in the resampling code is in c++ which goes through this switch statement to sort out how to map then numpy array -> c++ types

matplotlib/src/_image_wrapper.cpp

Lines 273 to 327 in ce334d4

	switch (PyArray_TYPE(input_array)) {
	case NPY_DOUBLE:
	Py_BEGIN_ALLOW_THREADS
	resample(
	(double *)PyArray_DATA(input_array),
	PyArray_DIM(input_array, 1),
	PyArray_DIM(input_array, 0),
	(double *)PyArray_DATA(output_array),
	PyArray_DIM(output_array, 1),
	PyArray_DIM(output_array, 0),
	params);
	Py_END_ALLOW_THREADS
	break;
	case NPY_FLOAT:
	Py_BEGIN_ALLOW_THREADS
	resample(
	(float *)PyArray_DATA(input_array),
	PyArray_DIM(input_array, 1),
	PyArray_DIM(input_array, 0),
	(float *)PyArray_DATA(output_array),
	PyArray_DIM(output_array, 1),
	PyArray_DIM(output_array, 0),
	params);
	Py_END_ALLOW_THREADS
	break;
	case NPY_UINT8:
	case NPY_BYTE:
	Py_BEGIN_ALLOW_THREADS
	resample(
	(unsigned char *)PyArray_DATA(input_array),
	PyArray_DIM(input_array, 1),
	PyArray_DIM(input_array, 0),
	(unsigned char *)PyArray_DATA(output_array),
	PyArray_DIM(output_array, 1),
	PyArray_DIM(output_array, 0),
	params);
	Py_END_ALLOW_THREADS
	break;
	case NPY_UINT16:
	case NPY_INT16:
	Py_BEGIN_ALLOW_THREADS
	resample(
	(unsigned short *)PyArray_DATA(input_array),
	PyArray_DIM(input_array, 1),
	PyArray_DIM(input_array, 0),
	(unsigned short *)PyArray_DATA(output_array),
	PyArray_DIM(output_array, 1),
	PyArray_DIM(output_array, 0),
	params);
	Py_END_ALLOW_THREADS
	break;
	default:
	PyErr_SetString(PyExc_ValueError, "Unsupported dtype");
	goto error;
	}

which does not support half-width floats (https://docs.scipy.org/doc/numpy/reference/c-api.dtype.html#c.NPY_HALF is not covered). The best fix for this would be to sort out what type we can use here. The underlying c++ code is all templated so if we cast the point to the correct type it should "just work", but I am not sure what that type is. I suspect we can get it from the numpy headers?

If we can't make that work, the python side fix is to up-cast the input data to np.float32.

[QuLogic]

I don't think there is a 16-bit floating-point type in C++, unless you use something like Half, so casting from NumPy's side is likely the only solution.

[QuLogic]

For reference, here's the information from NumPy: https://numpy.org/doc/1.18/reference/c-api/coremath.html#half-precision-functions

It's only a pseudo-type, and for any real math operations, the recommendation is to cast to float/double.