2

I have three numpy arrays, X, Y, and Z.

X and Y are coordinates of a spatial grid and each grid point (X, Y) has an intensity Z. I would like to save a PNG image using this data. Interpolation is not needed, as X and Y are guaranteed to cover each grid point between min(X) and max(Y).

I'm guessing the solution lies within numpy's meshgrid() function, but I can't figure out how to reshape the Z array to NxM intensity data.

How can I do that?

To clarify the input data structure, this is what it looks like:

   X   |    Y    |    Z
-----------------------------
0.1    | 0.1     | something..
0.1    | 0.2     | something..
0.1    | 0.3     | something..
...
0.2    | 0.1     | something..
0.2    | 0.2     | something..
0.2    | 0.3     | something..
...

0.2    | 0.1     | something..
0.1    | 0.2     | something..
0.3    | 0.3     | something..
...
MSeifert
  • 145,886
  • 38
  • 333
  • 352
janoliver
  • 7,744
  • 14
  • 60
  • 103

3 Answers3

2

To begin with, you should run this piece of code:

import numpy as np

X = np.asarray(<X data>)
Y = np.asarray(<Y data>)
Z = np.asarray(<Z data>)

Xu = np.unique(X)
Yu = np.unique(Y)

Then you could apply any of the following approaches. It is worth noting that all of them would work fine even if the data are NOT sorted (in contrast to the currently accepted answer):

1) A for loop and numpy.where() function

This is perhaps the simplest and most readable solution:

Zimg = np.zeros((Xu.size, Yu.size), np.uint8)
for i in range(X.size):
    Zimg[np.where(Xu==X[i]), np.where(Yu==Y[i])] = Z[i]

2) A list comprehension and numpy.sort() funtion

This solution - which is a bit more involved than the previous one - relies on Numpy's structured arrays:

data_type = [('x', np.float), ('y', np.float), ('z', np.uint8)]
XYZ = [(X[i], Y[i], Z[i]) for i in range(len(X))]
table = np.array(XYZ, dtype=data_type)
Zimg = np.sort(table, order=['y', 'x'])['z'].reshape(Xu.size, Yu.size)

3) Vectorization

Using lexsort is an elegant and efficient way of performing the required task:

Zimg = Z[np.lexsort((Y, X))].reshape(Xu.size, Yu.size)

4) Pure Python, not using NumPy

You may want to check out this link for a pure Python solution without any third party dependencies.

To end up, you have different options to save Zimg as an image:

from PIL import Image
Image.fromarray(Zimg).save('z-pil.png')

import matplotlib.pyplot as plt
plt.imsave('z-matplotlib.png', Zimg)

import cv2
cv2.imwrite('z-cv2.png', Zimg)

import scipy.misc
scipy.misc.imsave('z-scipy.png', Zimg)
Community
  • 1
  • 1
Tonechas
  • 13,398
  • 16
  • 46
  • 80
  • Hi, Z is a one dimensional array, just like X and Y. – janoliver Apr 29 '16 at 09:54
  • Thank you for your suggestion. I guess that would work, however, I was hoping for a more _pythonic_ way of doing it. But until I find something suitable, I'll stick to your code. – janoliver Apr 29 '16 at 10:36
  • 1
    If there is a large amount of data, consider using Pandas groupby method to speed up this operation: http://pandas.pydata.org/pandas-docs/stable/groupby.html – Peter Apr 29 '16 at 23:34
1

You said you needed no interpolation is needed since every grid point is covered. So I assume the points are equally spaced.

If your table is already sorted primary by increasing x and secondary by y you can simply take the Z array directly and save it using PIL:

import numpy as np
# Find out what shape your final array has (if you already know just hardcode these)
x_values = np.unique(X).size 
y_values = np.unique(Y).size 
img = np.reshape(Z, (x_values, y_values))

# Maybe you need to cast the dtype to fulfill png restrictions
#img = img.astype(np.uint) # alter it if needed

# Print image
from PIL import Image
Image.fromarray(img).save('filename.png')

If your input isn't sorted (it looks like it is but who knows) you have to sort it before you start. Depending on your input this can be easy or really hard.

MSeifert
  • 145,886
  • 38
  • 333
  • 352
  • reshape() is what I was looking for. Perfect, now it's a one-liner as I can simply do `np.unique(x).size`... Thank you! – janoliver Apr 29 '16 at 16:30
  • I can't understand how this code works. `y_values` is a tuple, how reshape can manage that ??? – B. M. Apr 29 '16 at 19:32
  • Thanks. But I think it works only if Z.size =x_values*y_values, what is not garanteed if there are duplicates..... – B. M. Apr 29 '16 at 19:48
  • From the question _"each grid point (X, Y) has an intensity Z"_ and _"X and Y are guaranteed to cover each grid point between min(X) and max(Y)"_ so it sounds like there are no duplicates and no missing values. – MSeifert Apr 29 '16 at 19:49
  • Ok, probably for the OP. that's not the case in the sample data. – B. M. Apr 29 '16 at 19:52
  • I didn't see that. Probably a copy mistake too. :-) – MSeifert Apr 29 '16 at 19:56
  • In the comments the author specifies that X and Y can have duplicate values. – Peter Apr 29 '16 at 23:33
  • @Peter I think he only meant that `X` for example contains multiple `0.1` and not that there are any `(X, Y)` duplicates. But that should not be questions for me, these questions should be for the OP. If there are duplicates and the values are unordered the answer needs to be different and if required I can also expand to give them but it already worked so why bother? To create a canonical answer? – MSeifert Apr 30 '16 at 11:40
0

np.ufunc.at is a good tool to manage duplicates in vectorized way.

Suppose these data :

In [3]: X,Y,Z=rand(3,10).round(1)

(array([ 0.4,  0.2,  0.1,  0.8,  0.4,  0.1,  0.5,  0.2,  0.6,  0.2]),
 array([ 0.5,  0.3,  0.5,  0.9,  0.9,  0.5,  0.3,  0.6,  0.4,  0.4]),
 array([ 0.4,  0.6,  0.6,  0.4,  0.1,  0.1,  0.2,  0.6,  0.9,  0.8]))

First scale the image (scale=3 here) :

In [4]: indices=[ (3*c).astype(int) for c in (X,Y)]

[array([1, 0, 0, 2, 1, 0, 1, 0, 1, 0]), array([1, 0, 1, 2, 2, 1, 0, 1, 1, 1])]

Make a empty image : image=zeros((3,3)), according to indices bounds.

Then build. Here we keep the maximum.

In [5]: np.maximum.at(image,indices,Z)  # in place

array([[ 0.6,  0.8,  0. ],
       [ 0.2,  0.9,  0.1],
       [ 0. ,  0. ,  0.4]])

Finally save in PNG : matplotlib.pyplot.imsave('img.png',image)

B. M.
  • 18,243
  • 2
  • 35
  • 54