SpyFile indexing result doesn't exactly match numpy indexing

[donm]

When indexing a numpy array, if a single index is specified along an axis (instead of a slice), then numpy removes the length-1 entries from the shape of the result. The __getitem__ function in SpyFile and some of interleave-specific routines (read_subregion, read_subimage) do not currently follow this convention.

Here are some examples that show the inconsistencies.

In [141]: spy.shape # spy is a SpyFile object (a BsqFile)
Out[141]: (640, 160, 6)

In [142]: arr = spy[:,:,:]; arr.shape # arr is a numpy array with the SpyFile data
Out[142]: (640, 160, 6)

In [143]: spy[:, :, 0].shape
Out[143]: (640, 160, 1)

In [144]: arr[:, :, 0].shape
Out[144]: (640, 160)

In [145]: arr[:, :, 0:1].shape                                                        
Out[145]: (640, 160, 1)

In [146]: spy.read_band(0).shape
Out[146]: (640, 160)

[tboggs]

The behavior in your examples is by design so I wouldn't consider it a bug. The web site documentation on the SpyFile Interface states that the SpyFile subscript (__getitem__) operator behaves "much like" numpy's and even shows an example like line 143 in your session above, though it would probably be better to have an example in the relevant doc strings. I'm guessing there are probably some examples of fancy indexing that work with numpy arrays but not SpyFile objects. So again, it is "much like" but not exactly the same.

I guess the real question here is whether not collapsing dimensions like numpy is a bad design choice and should be changed. My thinking when creating the __getitem__ interface for SpyFile objects was that people would primarily be using the interface to retrieve the data in an image-like layout. When reading a single band, it would often be convenient to remove the last dimension but the problem is when reading a single row or column of data.

Consider the following image:

In [26]: image = open_image('cup95eff.int.hdr')

In [27]: image.shape
Out[27]: (512, 614, 50)

In [28]: arr = image[:, :, :]   # or use image.load()

In [29]: arr.shape
Out[29]: (512, 614, 50)

You've shown above how the shapes work out when reading a single band. Now suppose that I want to pass a column of the image to a function that operates on an image. Consider the shapes of the column image data for the two cases:

In [30]: image[:, 10, :].shape
Out[30]: (512, 1, 50)

In [31]: arr[:, 10, :].shape
Out[31]: (512, 50)

The SpyFile interface gives an array of image data with 512 rows, 1 column, and 50 bands. But what does the shape (512, 50) mean. A single-row image, a single-column image, or a single-band image? Usually, a shape (M, N) would be interpreted by a function as a single-band image but that isn't the case here. That's basically why I made the interface the way I did.

When reading a single band (i.e., image[:, :, b]), it would usually be convenient to collapse the last dimension (like the numpy array does) but then I would have different rules for how the last dimension is handled relative to the first two and I thought that would be overly complicated. So that is why the read_band method exists (I did actually remember to note that in the doc string).

This all said, if you think this feature of the interface makes it overly complicated or error prone, I'm open to weighing the options of changing it.

[donm]

I was wondering if it was intentional and these reasons make sense. I had guessed some of them but my gut reaction was still, "but...consistency..." But I hadn't thought about how it is easier to implement a function if you can always expect to have a cube of data instead of maybe a cube, or a slice that isn't a band, or a pixel. I think that that alone probably makes it worth the convention of always dealing with cubes even when slicing.

Plus, I guess it's not that big of a deal to just add .squeeze() when you want a 2D slice.

[tboggs]

I suppose there could be a settable option in the SpyFile class to exactly mimic numpy's indexing behavior. That is fairly easy to do when the SpyFile object uses a memmap object to access the data because we could do something like the following in the SpyFile subclass __getitem__ method:

data = self._memmap.__getitem__(*kwargs)

(though there would be the additional details of handling data interleave formats)

The harder part would be in cases where the SpyFile isn't using a memmap (e.g., when accessing a very large image file on a memory constrained system). SpyFile.__getitem__ converts it's arguments into lists of index values that are passed to read_subimage and that method wouldn't know if its third argument (spectral dimension) value of [3] was the result of SpyFile.__getitem__ receiving 3 or slice(3, 4) as the corresponding value.

But the more I think about this, I'm inclined to agree that it's not a big deal to add .squeeze().

I'll leave this issue open for a while in case anyone has any additional ideas or considerations.

[donm]

We've been talking about SpyFile vs numpy arrays, but do you think that it's a problem that SpyFile slicing has different behavior compared to an ImageArray? If the SpyFile behavior stays the same, then I think (without having thought it through all the way) I might advocate for changing the ImageArray slicing behavior to match.

In [98]: spy.shape                                                                                                                                                                                                  
Out[98]: (500, 128, 128)

In [99]: arr = spy.load()

In [100]: spy[:,:,0].shape                                                                                                                                                                                          
Out[100]: (500, 128, 1) # if this is how SpyFile behaves...

In [101]: arr[:,:,0].shape                                                                                                                                                                                          
Out[101]: (500, 128)  # ...then I think I would consider this unexpected

What do you think? This change would probably be easy to make by intercepting the ImageArray.__getitem__ call and replacing any single indices with a slice (unless they're all single) and then passing the arguments to the parent __getitem__.

[tboggs]

Good point. That is an internal inconsistency and should be fixed. I think a user should reasonable expect an ImageArray to behave like a SpyFile object in that sense. Overloading ImageArray.__getitem__ as you mentioned should work.

[donm]

I'm working on this, by the way, and will have something to you tomorrow. Or rather, later today.