backwards incompatibility: algebraic operations with numpy arrays/scalars on the "left" return a numpy array

[KrisThielemans]

In the previous version of SIRF (e.g. used for PETRIC1), things like sirf_image * numpy_array worked, but now it raises an error. This breaks for instance the MaGeZ implementation. https://github.com/SyneRBI/PETRIC2/issues/12.

@evgueni-ovtchinnikov is it possible to restore this?

[evgueni-ovtchinnikov]

I firmly believe that allowing things like a*b where a is a SIRF object and b is neither a SIRF object nor a scalar is not a good idea, as it is likely to do "papering over the crack" (the actual mistake of using same-named variable for two objects of different type was done elsewhere and must be corrected).

[KrisThielemans]

@evgueni-ovtchinnikov Actually, the underlying (and main) problem is that currently number * STIR.ImageData returns a numpy array. Not sure if #1358 fixes that? would be best to add a test for it.

[evgueni-ovtchinnikov]

There appears to be a STIR issue: I put this debugging print in stir_data_containers.h:

		virtual bool supports_array_view() const
		{
#if STIR_VERSION >= 060200
			std::cout << data().is_contiguous() << "!!!\n";
			return data().is_contiguous();
#else
			return false;
#endif
		}

and when running this in a Python script:

    image = acq_data.create_uniform_image(1.0)
    img_arr = image.as_array()
    x = img_arr*image
    print(type(x))

I got

0!!!

To remind you:

		Image3DF& data()
		{
			return *_data;
		}
		const Image3DF& data() const
		{
			return *_data;
		}

and

		stir::shared_ptr<Image3DF> _data;

and

	typedef stir::DiscretisedDensity<3, float> Image3DF;

[KrisThielemans]

not sure which of all of these objects is printing the 0!!!. Can you split this up a bit? Also, I'm assuming SIRF is handling the numerical operations first (although possibly by just forwarding to STIR).

[KrisThielemans]

The STIR issue turns out to be independent of the main issue here.

We believe that numpy array algebraic operations such as __mul__ internally call asarray() on the rhs argument, which now works since SIRF 3.9. Due to Python rules, this means that we cannot "intercept" this call. We therefore see no way to fix this in SIRF, and will address this by documenting this behaviour in UsersGuide.md.

[paskino]

I think our unit tests do not test for inplace algebric operation as +=. These should be added in the DataContainerAlgebraTests class and then they'll run for all data container type.