perf: oindex optimization (zarr-developers#3830)

* oindex single dim optimization

* changelog

* changelog type

* Apply suggestions from code review

Co-authored-by: Deepak Cherian <dcherian@users.noreply.github.com>

* lint

---------

Co-authored-by: Deepak Cherian <dcherian@users.noreply.github.com>

Author: maxrjones

changes/3830.misc.md

@@ -0,0 +1 @@
+Optimize the performance of indexing operations when using an array-like indexer on a single dimension.

src/zarr/core/indexing.py

@@ -981,19 +981,26 @@ def __iter__(self) -> Iterator[ChunkProjection]:
 
             # handle advanced indexing arrays orthogonally
             if self.is_advanced:
-                # N.B., numpy doesn't support orthogonal indexing directly as yet,
-                # so need to work around via np.ix_. Also np.ix_ does not support a
-                # mixture of arrays and slices or integers, so need to convert slices
-                # and integers into ranges.
-                chunk_shape = tuple(
-                    g.chunk_size(p.dim_chunk_ix)
-                    for g, p in zip(self.dim_grids, dim_projections, strict=True)
-                )
-                chunk_selection = ix_(chunk_selection, chunk_shape)
+                # NumPy can handle a single array-indexed dimension directly,
+                # which preserves full slices and avoids an
+                # unnecessary advanced-indexing copy. Integer-indexed
+                # dimensions still need the ix_ path for downstream squeezing.
+                # Example: we skip `ix_` for array[:, :, [1, 2, 3]]
+                n_array_dims = sum(isinstance(sel, np.ndarray) for sel in chunk_selection)
+
+                if n_array_dims > 1 or self.drop_axes:
+                    # N.B., numpy doesn't support orthogonal indexing directly
+                    # for multiple array-indexed dimensions, so we need to
+                    # convert the orthogonal selection into coordinate arrays.
+                    chunk_shape = tuple(
+                        g.chunk_size(p.dim_chunk_ix)
+                        for g, p in zip(self.dim_grids, dim_projections, strict=True)
+                    )
+                    chunk_selection = ix_(chunk_selection, chunk_shape)
 
-                # special case for non-monotonic indices
-                if not is_basic_selection(out_selection):
-                    out_selection = ix_(out_selection, self.shape)
+                    # special case for non-monotonic indices
+                    if not is_basic_selection(out_selection):
+                        out_selection = ix_(out_selection, self.shape)
 
             is_complete_chunk = all(p.is_complete_chunk for p in dim_projections)
             yield ChunkProjection(chunk_coords, chunk_selection, out_selection, is_complete_chunk)