1

I have a 3D numpy array that I want to rotate with an angle that I want. I have tried using scipy.ndimage.rotate function and it does the job. However, it does a lot of rounding when rotating. This causes me a problem because my 3D array is representation of an object and numbers in each pixel represent the material that pixel is filled with (which I store in a different file). Therefore, I need a way to rotate the array without doing approximation or rounding and making the object blurry is not a problem

Here is what I got with the function I used:

Before

After

Toykan Ozdeger
  • 139
  • 2
  • 3
  • 12
  • 1
    rotating around 90 degrees? Can you give an example of before and after the tranformation? – some_name.py Jun 19 '19 at 13:13
  • @user7784503 No 90 degrees is easy. The problem is the angles between. I have updated the question to show the results – Toykan Ozdeger Jun 19 '19 at 13:19
  • I think the key is your interpolation method, you want to find an implementation that uses nearest neighbor interpolation rather than spline / cubic etc. Is that correct? – brechmos Jun 19 '19 at 13:20
  • You might want to look at https://codereview.stackexchange.com/questions/41688/rotating-greyscale-images (well, that is 2D, but you could maybe extend the idea) – brechmos Jun 19 '19 at 13:22
  • @brechmos Yes that's correct – Toykan Ozdeger Jun 19 '19 at 13:25
  • @brechmos The code in the link is incomplete I think. there is a link below the code that says "full code" but the link doesn't work. – Toykan Ozdeger Jun 19 '19 at 13:54

3 Answers3

2

The problem you are dealing with is essentially a sampling issue. Your resolution is too low for the data you are dealing with. One possibility to solve this is to increase the resolution of the image you are working with, enforce the color values as you rotate (ie no blending colors at the edges), and create a size/shape template that must be met after the rotation.

Edit: For clarity, it isn't the data that is at too low of a resolution, it's the image in which the data is stored that should be at a high enough resolution. The wikipedia page on multidimensional sampling is good for this topic: https://en.wikipedia.org/wiki/Multidimensional_sampling

  • The image is result of a simulation. This is why the resolution is low and it isn't possible to get a better resolution. I exactly need to do what you explained I think but just no idea how to. – Toykan Ozdeger Jun 19 '19 at 13:58
  • Without knowing exactly how you're doing this I can't say exactly how to do it, but this question has good answers I think: https://stackoverflow.com/questions/48121916/numpy-resize-rescale-image –  Jun 19 '19 at 14:20
1

I think the way I would approach it, outside of someone knowing an actual package to do this, is start with the indices and rotate them, then, given they may be floating point, round them. This may not be the best, but I think it should work.

Most of this example is loading a 3D dataset I found to use as an example.

import matplotlib.pyplot as plt
import os
import numpy as np
from scipy.ndimage import rotate

def load_example_data():
    # Found data as an example
    from urllib.request import urlopen
    import tarfile
    opener = urlopen( 'http://graphics.stanford.edu/data/voldata/MRbrain.tar.gz')
    tar_file = tarfile.open('MRbrain.tar.gz')
    try:
        os.mkdir('mri_data')
    except:
        pass
    tar_file.extractall('mri_data')
    tar_file.close()

    import numpy as np
    data = np.array([np.fromfile(os.path.join('mri_data', 'MRbrain.%i' % i),
                                            dtype='>u2') for i in range(1, 110)])
    data.shape = (109, 256, 256)

    return data

def rotate_nn(data, angle, axes):
    """
    Rotate a `data` based on rotating coordinates.
    """

    # Create grid of indices
    shape = data.shape
    d1, d2, d3 = np.mgrid[0:shape[0], 0:shape[1], 0:shape[2]]

    # Rotate the indices
    d1r = rotate(d1, angle=angle, axes=axes)
    d2r = rotate(d2, angle=angle, axes=axes)
    d3r = rotate(d3, angle=angle, axes=axes)

    # Round to integer indices
    d1r = np.round(d1r)
    d2r = np.round(d2r)
    d3r = np.round(d3r)

    d1r = np.clip(d1r, 0, shape[0])
    d2r = np.clip(d2r, 0, shape[1])
    d3r = np.clip(d3r, 0, shape[2])

    return data[d1r, d2r, d3r]


data = load_example_data()

# Rotate the coordinates indices
angle = 5
axes = (0, 1)

data_r = rotate_nn(data, angle, axes)

I think the general idea will work. You will have to consider what the axis is to rotate around.

brechmos
  • 1,278
  • 1
  • 11
  • 22
  • Hmmm... come to think of it, you might have to consider indices that get rotated out of the `0..data.shape` range. So, might need to add an `np.clip()` after `rotate` and`np.round` in the `rotate_nn` method. I'll see if I can clean that up. – brechmos Jun 19 '19 at 14:26
  • ValueError: The truth value of an array with more than one element is ambiguous. Use a.any() or a.all(). I seem to get this error with the above code – Toykan Ozdeger Jun 19 '19 at 16:21
  • Is the error when you use your data, or the example data above? – brechmos Jun 19 '19 at 16:22
  • yes it is when I use my data. it is a 60,60,25 array with integers only – Toykan Ozdeger Jun 19 '19 at 16:24
  • Hmmm... I would suggest you check the d1 and compare it to d1r (same for d2/d3) and make sure the rotation isn't doing something odd. (Might be the `axes` parameter?). Not sure without working with the actual data. – brechmos Jun 19 '19 at 16:33
  • d1r, d2r and d3r are all 3D arrays I think. And the error is because we are trying to put three 3D arrays into one array – Toykan Ozdeger Jun 19 '19 at 17:37
  • Let us [continue this discussion in chat](https://chat.stackoverflow.com/rooms/195229/discussion-between-brechmos-and-toykan-ozdeger). – brechmos Jun 20 '19 at 00:24
1

For anyone with this problem stumbling upon this thread: brechmos' comment under the OP put me in the right direction for an actual solution. rotate() by default uses a third-order spline interpolation, which gives nice smooth edges. We want sharp edges though, without numbers in between. Setting order = 0 does exactly this. No need for extra functions or implementing anything yourself, just change a single argument.

Bartian
  • 33
  • 5