ENH: missingdata: Make arr.item() and arr.itemset() work with NA masks

5 files changed, 408 insertions, 193 deletions

diff --git a/numpy/core/arrayprint.py b/numpy/core/arrayprint.py
index 1d072fe5d..506cce8cc 100644
--- a/numpy/core/arrayprint.py
+++ b/numpy/core/arrayprint.py
@@ -430,16 +430,20 @@ def array2string(a, max_line_width=None, precision=None,
 
     if a.shape == ():
         x = a.item()
-        try:
-            lst = a._format(x)
-            msg = "The `_format` attribute is deprecated in Numpy 2.0 and " \
-                  "will be removed in 2.1. Use the `formatter` kw instead."
-            import warnings
-            warnings.warn(msg, DeprecationWarning)
-        except AttributeError:
-            if isinstance(x, tuple):
-                x = _convert_arrays(x)
-            lst = style(x)
+        if isna(x):
+            lst = str(x)
+        else:
+            try:
+                lst = a._format(x)
+                msg = "The `_format` attribute is deprecated in Numpy " \
+                      "2.0 and will be removed in 2.1. Use the " \
+                      "`formatter` kw instead."
+                import warnings
+                warnings.warn(msg, DeprecationWarning)
+            except AttributeError:
+                if isinstance(x, tuple):
+                    x = _convert_arrays(x)
+                lst = style(x)
     elif reduce(product, a.shape) == 0:
         # treat as a null array if any of shape elements == 0
         lst = "[]"
diff --git a/numpy/core/src/multiarray/mapping.c b/numpy/core/src/multiarray/mapping.c
index dc1dd66f4..f016cd755 100644
--- a/numpy/core/src/multiarray/mapping.c
+++ b/numpy/core/src/multiarray/mapping.c
@@ -466,7 +466,7 @@ count_new_axes_0d(PyObject *tuple)
 NPY_NO_EXPORT PyObject *
 add_new_axes_0d(PyArrayObject *arr,  int newaxis_count)
 {
-    PyArrayObject *other;
+    PyArrayObject *ret;
     npy_intp dimensions[NPY_MAXDIMS];
     int i;
 
@@ -474,21 +474,41 @@ add_new_axes_0d(PyArrayObject *arr,  int newaxis_count)
         dimensions[i]  = 1;
     }
     Py_INCREF(PyArray_DESCR(arr));
-    other = (PyArrayObject *)PyArray_NewFromDescr(Py_TYPE(arr),
+    ret = (PyArrayObject *)PyArray_NewFromDescr(Py_TYPE(arr),
                                 PyArray_DESCR(arr),
                                 newaxis_count, dimensions,
                                 NULL, PyArray_DATA(arr),
-                                PyArray_FLAGS(arr),
+            PyArray_FLAGS(arr) & ~(NPY_ARRAY_MASKNA | NPY_ARRAY_OWNMASKNA),
                                 (PyObject *)arr);
-    if (other == NULL) {
+    if (ret == NULL) {
         return NULL;
     }
+
     Py_INCREF(arr);
-    if (PyArray_SetBaseObject(other, (PyObject *)arr) < 0) {
-        Py_DECREF(other);
+    if (PyArray_SetBaseObject(ret, (PyObject *)arr) < 0) {
+        Py_DECREF(ret);
         return NULL;
     }
-    return (PyObject *)other;
+
+    /* Take a view of the NA mask if it exists */
+    if (PyArray_HASMASKNA(arr)) {
+        PyArrayObject_fieldaccess *fret = (PyArrayObject_fieldaccess *)ret;
+
+        fret->maskna_dtype = PyArray_MASKNA_DTYPE(arr);
+        Py_INCREF(fret->maskna_dtype);
+
+        fret->maskna_data = PyArray_MASKNA_DATA(arr);
+
+        for (i = 0; i < newaxis_count; ++i) {
+            fret->maskna_strides[i]  = fret->maskna_dtype->elsize;
+        }
+
+        /* This view doesn't own the mask */
+        fret->flags |= NPY_ARRAY_MASKNA;
+        fret->flags &= ~NPY_ARRAY_OWNMASKNA;
+    }
+
+    return (PyObject *)ret;
 }
 
 
@@ -1297,7 +1317,8 @@ array_subscript(PyArrayObject *self, PyObject *op)
                 Py_INCREF(self);
                 return (PyObject *)self;
             }
-            if ((nd = count_new_axes_0d(op)) == -1) {
+            nd = count_new_axes_0d(op);
+            if (nd == -1) {
                 return NULL;
             }
             return add_new_axes_0d(self, nd);
@@ -1776,18 +1797,20 @@ array_subscript_nice(PyArrayObject *self, PyObject *op)
      * array_subscript_simple).  So, this cast is a bit dangerous..
      */
 
-    /*
-     * The following is just a copy of PyArray_Return with an
-     * additional logic in the nd == 0 case.
-     */
-
     if (mp == NULL) {
         return NULL;
     }
+
     if (PyErr_Occurred()) {
         Py_XDECREF(mp);
         return NULL;
     }
+
+    /*
+     * The following adds some additional logic to avoid calling
+     * PyArray_Return if there is an ellipsis.
+     */
+
     if (PyArray_Check(mp) && PyArray_NDIM(mp) == 0) {
         npy_bool noellipses = TRUE;
         if ((op == Py_Ellipsis) || PyString_Check(op) || PyUnicode_Check(op)) {
@@ -1815,12 +1838,10 @@ array_subscript_nice(PyArrayObject *self, PyObject *op)
             }
         }
         if (noellipses) {
-            PyObject *ret;
-            ret = PyArray_ToScalar(PyArray_DATA(mp), mp);
-            Py_DECREF(mp);
-            return ret;
+            return PyArray_Return(mp);
         }
     }
+
     return (PyObject *)mp;
 }
 
diff --git a/numpy/core/src/multiarray/methods.c b/numpy/core/src/multiarray/methods.c
index c728f17e0..3446ec738 100644
--- a/numpy/core/src/multiarray/methods.c
+++ b/numpy/core/src/multiarray/methods.c
@@ -19,6 +19,7 @@
 
 #include "methods.h"
 #include "convert_datatype.h"
+#include "na_singleton.h"
 
 
 /* NpyArg_ParseKeywords
@@ -422,7 +423,7 @@ NPY_NO_EXPORT PyObject *
 PyArray_Byteswap(PyArrayObject *self, Bool inplace)
 {
     PyArrayObject *ret;
-    intp size;
+    npy_intp size;
     PyArray_CopySwapNFunc *copyswapn;
     PyArrayIterObject *it;
 
@@ -440,7 +441,7 @@ PyArray_Byteswap(PyArrayObject *self, Bool inplace)
         }
         else { /* Use iterator */
             int axis = -1;
-            intp stride;
+            npy_intp stride;
             it = (PyArrayIterObject *)                      \
                 PyArray_IterAllButAxis((PyObject *)self, &axis);
             stride = PyArray_STRIDES(self)[axis];
@@ -549,103 +550,251 @@ array_tofile(PyArrayObject *self, PyObject *args, PyObject *kwds)
     return Py_None;
 }
 
-
+/*
+ * Gets a single item from the array, based on a single multi-index
+ * array of values, which must be of length PyArray_NDIM(self).
+ */
 static PyObject *
-array_toscalar(PyArrayObject *self, PyObject *args) {
-    int n, nd;
-    n = PyTuple_GET_SIZE(args);
+PyArray_MultiIndexGetItem(PyArrayObject *self, npy_intp *multi_index)
+{
+    int idim, ndim = PyArray_NDIM(self);
+    char *data = PyArray_DATA(self);
+    npy_intp *shape = PyArray_SHAPE(self);
+    npy_intp *strides = PyArray_STRIDES(self);
 
-    if (n == 1) {
-        PyObject *obj;
-        obj = PyTuple_GET_ITEM(args, 0);
-        if (PyTuple_Check(obj)) {
-            args = obj;
-            n = PyTuple_GET_SIZE(args);
+    /* Case with an NA mask */
+    if (PyArray_HASMASKNA(self)) {
+        char *maskdata = PyArray_MASKNA_DATA(self);
+        npy_mask maskvalue;
+        npy_intp *maskstrides = PyArray_MASKNA_STRIDES(self);
+
+        if (PyArray_HASFIELDS(self)) {
+            PyErr_SetString(PyExc_RuntimeError,
+                    "field-NA is not supported yet in MultiIndexGetItem");
+            return NULL;
         }
-    }
 
-    if (n == 0) {
-        if (PyArray_NDIM(self) == 0 || PyArray_SIZE(self) == 1)
-            return PyArray_DESCR(self)->f->getitem(PyArray_DATA(self), self);
+        /* Get the data and maskdata pointer */
+        for (idim = 0; idim < ndim; ++idim) {
+            npy_intp shapevalue = shape[idim];
+            npy_intp ind = multi_index[idim];
+
+            if (ind < 0) {
+                ind += shapevalue;
+            }
+
+            if (ind < 0 || ind >= shapevalue) {
+                PyErr_SetString(PyExc_ValueError, "index out of bounds");
+                return NULL;
+            }
+
+            data += ind * strides[idim];
+            maskdata += ind * maskstrides[idim];
+        }
+
+        maskvalue = (npy_mask)*maskdata;
+        if (NpyMaskValue_IsExposed(maskvalue)) {
+            return PyArray_DESCR(self)->f->getitem(data, self);
+        }
         else {
-            PyErr_SetString(PyExc_ValueError,
-                            "can only convert an array "    \
-                            " of size 1 to a Python scalar");
-            return NULL;
+            return (PyObject *)NpyNA_FromDTypeAndMaskValue(
+                                                PyArray_DTYPE(self),
+                                                maskvalue, 0);
         }
     }
-    else if (n != PyArray_NDIM(self) && (n > 1 || PyArray_NDIM(self) == 0)) {
-        PyErr_SetString(PyExc_ValueError,
-                        "incorrect number of indices for "      \
-                        "array");
-        return NULL;
+    /* Case without an NA mask */
+    else {
+        /* Get the data pointer */
+        for (idim = 0; idim < ndim; ++idim) {
+            npy_intp shapevalue = shape[idim];
+            npy_intp ind = multi_index[idim];
+
+            if (ind < 0) {
+                ind += shapevalue;
+            }
+
+            if (ind < 0 || ind >= shapevalue) {
+                PyErr_SetString(PyExc_ValueError, "index out of bounds");
+                return NULL;
+            }
+
+            data += ind * strides[idim];
+        }
+
+        return PyArray_DESCR(self)->f->getitem(data, self);
     }
-    else if (n == 1) { /* allows for flat getting as well as 1-d case */
-        intp value, loc, index, factor;
-        intp factors[MAX_DIMS];
-        value = PyArray_PyIntAsIntp(PyTuple_GET_ITEM(args, 0));
-        if (error_converting(value)) {
-            PyErr_SetString(PyExc_ValueError, "invalid integer");
-            return NULL;
+}
+
+/*
+ * Sets a single item in the array, based on a single multi-index
+ * array of values, which must be of length PyArray_NDIM(self).
+ *
+ * Returns 0 on success, -1 on failure.
+ */
+static int
+PyArray_MultiIndexSetItem(PyArrayObject *self, npy_intp *multi_index,
+                                                PyObject *obj)
+{
+    int idim, ndim = PyArray_NDIM(self);
+    char *data = PyArray_DATA(self);
+    npy_intp *shape = PyArray_SHAPE(self);
+    npy_intp *strides = PyArray_STRIDES(self);
+
+    /* Case with an NA mask */
+    if (PyArray_HASMASKNA(self)) {
+        char *maskdata = PyArray_MASKNA_DATA(self);
+        npy_intp *maskstrides = PyArray_MASKNA_STRIDES(self);
+        NpyNA *na = NpyNA_FromObject(obj, 1);
+
+        if (PyArray_HASFIELDS(self)) {
+            PyErr_SetString(PyExc_RuntimeError,
+                    "field-NA is not supported yet in MultiIndexSetItem");
+            return -1;
         }
-        factor = PyArray_SIZE(self);
-        if (value < 0) value += factor;
-        if ((value >= factor) || (value < 0)) {
-            PyErr_SetString(PyExc_ValueError,
-                            "index out of bounds");
-            return NULL;
+
+        /* Get the data and maskdata pointer */
+        for (idim = 0; idim < ndim; ++idim) {
+            npy_intp shapevalue = shape[idim];
+            npy_intp ind = multi_index[idim];
+
+            if (ind < 0) {
+                ind += shapevalue;
+            }
+
+            if (ind < 0 || ind >= shapevalue) {
+                PyErr_SetString(PyExc_ValueError, "index out of bounds");
+                return -1;
+            }
+
+            data += ind * strides[idim];
+            maskdata += ind * maskstrides[idim];
         }
-        if (PyArray_NDIM(self) == 1) {
-            value *= PyArray_STRIDES(self)[0];
-            return PyArray_DESCR(self)->f->getitem(PyArray_DATA(self) + value,
-                                           self);
+
+        if (na == NULL) {
+            *maskdata = 1;
+            return PyArray_DESCR(self)->f->setitem(obj, data, self);
         }
-        nd = PyArray_NDIM(self);
-        factor = 1;
-        while (nd--) {
-            factors[nd] = factor;
-            factor *= PyArray_DIMS(self)[nd];
+        else {
+            char maskvalue = (char)NpyNA_AsMaskValue(na);
+
+            if (maskvalue != 0 &&
+                        PyArray_MASKNA_DTYPE(self)->type_num != NPY_MASK) {
+                /* TODO: also handle struct-NA mask dtypes */
+                PyErr_SetString(PyExc_ValueError,
+                        "Cannot assign an NA with a payload to an "
+                        "NA-array with a boolean mask, requires a "
+                        "multi-NA mask");
+                return -1;
+            }
+
+            *maskdata = maskvalue;
+
+            return 0;
         }
-        loc = 0;
-        for (nd = 0; nd < PyArray_NDIM(self); nd++) {
-            index = value / factors[nd];
-            value = value % factors[nd];
-            loc += PyArray_STRIDES(self)[nd]*index;
+    }
+    /* Case without an NA mask */
+    else {
+        /* Get the data pointer */
+        for (idim = 0; idim < ndim; ++idim) {
+            npy_intp shapevalue = shape[idim];
+            npy_intp ind = multi_index[idim];
+
+            if (ind < 0) {
+                ind += shapevalue;
+            }
+
+            if (ind < 0 || ind >= shapevalue) {
+                PyErr_SetString(PyExc_ValueError, "index out of bounds");
+                return -1;
+            }
+
+            data += ind * strides[idim];
         }
 
-        return PyArray_DESCR(self)->f->getitem(PyArray_DATA(self) + loc,
-                                       self);
+        return PyArray_DESCR(self)->f->setitem(obj, data, self);
+    }
+}
+
 
+static PyObject *
+array_toscalar(PyArrayObject *self, PyObject *args)
+{
+    npy_intp multi_index[NPY_MAXDIMS];
+    int n = PyTuple_GET_SIZE(args);
+    int idim, ndim = PyArray_NDIM(self);
+
+    /* If there is a tuple as a single argument, treat it as the argument */
+    if (n == 1 && PyTuple_Check(PyTuple_GET_ITEM(args, 0))) {
+        args = PyTuple_GET_ITEM(args, 0);
+        n = PyTuple_GET_SIZE(args);
     }
-    else {
-        intp loc, index[MAX_DIMS];
-        nd = PyArray_IntpFromSequence(args, index, MAX_DIMS);
-        if (nd < n) {
+
+    if (n == 0) {
+        if (PyArray_SIZE(self) == 1) {
+            for (idim = 0; idim < ndim; ++idim) {
+                multi_index[idim] = 0;
+            }
+        }
+        else {
+            PyErr_SetString(PyExc_ValueError,
+                            "can only convert an array "
+                            " of size 1 to a Python scalar");
+        }
+    }
+    /* Special case of C-order flat indexing... :| */
+    else if (n == 1 && ndim != 1) {
+        npy_intp *shape = PyArray_SHAPE(self);
+        npy_intp value, size = PyArray_SIZE(self);
+
+        value = PyArray_PyIntAsIntp(PyTuple_GET_ITEM(args, 0));
+        if (value == -1 && PyErr_Occurred()) {
             return NULL;
         }
-        loc = 0;
-        while (nd--) {
-            if (index[nd] < 0) {
-                index[nd] += PyArray_DIMS(self)[nd];
-            }
-            if (index[nd] < 0 ||
-                index[nd] >= PyArray_DIMS(self)[nd]) {
-                PyErr_SetString(PyExc_ValueError,
-                                "index out of bounds");
+
+        /* Negative indexing */
+        if (value < 0) {
+            value += size;
+        }
+
+        if (value < 0 || value >= size) {
+            PyErr_SetString(PyExc_ValueError, "index out of bounds");
+            return NULL;
+        }
+
+        /* Convert the flat index into a multi-index */
+        for (idim = ndim-1; idim >= 0; --idim) {
+            multi_index[idim] = value % shape[idim];
+            value /= shape[idim];
+        }
+    }
+    /* A multi-index tuple */
+    else if (n == ndim) {
+        npy_intp value;
+
+        for (idim = 0; idim < ndim; ++idim) {
+            value = PyArray_PyIntAsIntp(PyTuple_GET_ITEM(args, idim));
+            if (value == -1 && PyErr_Occurred()) {
                 return NULL;
             }
-            loc += PyArray_STRIDES(self)[nd]*index[nd];
+            multi_index[idim] = value;
         }
-        return PyArray_DESCR(self)->f->getitem(PyArray_DATA(self) + loc, self);
     }
+    else {
+        PyErr_SetString(PyExc_ValueError,
+                        "incorrect number of indices for array");
+        return NULL;
+    }
+
+    return PyArray_MultiIndexGetItem(self, multi_index);
 }
 
 static PyObject *
-array_setscalar(PyArrayObject *self, PyObject *args) {
-    int n, nd;
-    int ret = -1;
+array_setscalar(PyArrayObject *self, PyObject *args)
+{
+    npy_intp multi_index[NPY_MAXDIMS];
+    int n = PyTuple_GET_SIZE(args) - 1;
+    int idim, ndim = PyArray_NDIM(self);
     PyObject *obj;
-    n = PyTuple_GET_SIZE(args) - 1;
 
     if (n < 0) {
         PyErr_SetString(PyExc_ValueError,
@@ -653,110 +802,76 @@ array_setscalar(PyArrayObject *self, PyObject *args) {
         return NULL;
     }
     obj = PyTuple_GET_ITEM(args, n);
+
+    /* If there is a tuple as a single argument, treat it as the argument */
+    if (n == 1 && PyTuple_Check(PyTuple_GET_ITEM(args, 0))) {
+        args = PyTuple_GET_ITEM(args, 0);
+        n = PyTuple_GET_SIZE(args);
+    }
+
     if (n == 0) {
-        if (PyArray_NDIM(self) == 0 || PyArray_SIZE(self) == 1) {
-            ret = PyArray_DESCR(self)->f->setitem(obj, PyArray_DATA(self), self);
+        if (PyArray_SIZE(self) == 1) {
+            for (idim = 0; idim < ndim; ++idim) {
+                multi_index[idim] = 0;
+            }
         }
         else {
             PyErr_SetString(PyExc_ValueError,
-                            "can only place a scalar for an "
-                            " array of size 1");
-            return NULL;
+                            "can only convert an array "
+                            " of size 1 to a Python scalar");
         }
     }
-    else if (n != PyArray_NDIM(self) && (n > 1 || PyArray_NDIM(self) == 0)) {
-        PyErr_SetString(PyExc_ValueError,
-                        "incorrect number of indices for "      \
-                        "array");
-        return NULL;
-    }
-    else if (n == 1) { /* allows for flat setting as well as 1-d case */
-        intp value, loc, index, factor;
-        intp factors[MAX_DIMS];
-        PyObject *indobj;
+    /* Special case of C-order flat indexing... :| */
+    else if (n == 1 && ndim != 1) {
+        npy_intp *shape = PyArray_SHAPE(self);
+        npy_intp value, size = PyArray_SIZE(self);
 
-        indobj = PyTuple_GET_ITEM(args, 0);
-        if (PyTuple_Check(indobj)) {
-            PyObject *res;
-            PyObject *newargs;
-            PyObject *tmp;
-            int i, nn;
-            nn = PyTuple_GET_SIZE(indobj);
-            newargs = PyTuple_New(nn+1);
-            Py_INCREF(obj);
-            for (i = 0; i < nn; i++) {
-                tmp = PyTuple_GET_ITEM(indobj, i);
-                Py_INCREF(tmp);
-                PyTuple_SET_ITEM(newargs, i, tmp);
-            }
-            PyTuple_SET_ITEM(newargs, nn, obj);
-            /* Call with a converted set of arguments */
-            res = array_setscalar(self, newargs);
-            Py_DECREF(newargs);
-            return res;
-        }
-        value = PyArray_PyIntAsIntp(indobj);
-        if (error_converting(value)) {
-            PyErr_SetString(PyExc_ValueError, "invalid integer");
-            return NULL;
-        }
-        if (value >= PyArray_SIZE(self)) {
-            PyErr_SetString(PyExc_ValueError,
-                            "index out of bounds");
+        value = PyArray_PyIntAsIntp(PyTuple_GET_ITEM(args, 0));
+        if (value == -1 && PyErr_Occurred()) {
             return NULL;
         }
-        if (PyArray_NDIM(self) == 1) {
-            value *= PyArray_STRIDES(self)[0];
-            ret = PyArray_DESCR(self)->f->setitem(obj, PyArray_DATA(self) + value,
-                                          self);
-            goto finish;
-        }
-        nd = PyArray_NDIM(self);
-        factor = 1;
-        while (nd--) {
-            factors[nd] = factor;
-            factor *= PyArray_DIMS(self)[nd];
-        }
-        loc = 0;
-        for (nd = 0; nd < PyArray_NDIM(self); nd++) {
-            index = value / factors[nd];
-            value = value % factors[nd];
-            loc += PyArray_STRIDES(self)[nd]*index;
+
+        /* Negative indexing */
+        if (value < 0) {
+            value += size;
         }
 
-        ret = PyArray_DESCR(self)->f->setitem(obj, PyArray_DATA(self) + loc, self);
-    }
-    else {
-        intp loc, index[MAX_DIMS];
-        PyObject *tupargs;
-        tupargs = PyTuple_GetSlice(args, 0, n);
-        nd = PyArray_IntpFromSequence(tupargs, index, MAX_DIMS);
-        Py_DECREF(tupargs);
-        if (nd < n) {
+        if (value < 0 || value >= size) {
+            PyErr_SetString(PyExc_ValueError, "index out of bounds");
             return NULL;
         }
-        loc = 0;
-        while (nd--) {
-            if (index[nd] < 0) {
-                index[nd] += PyArray_DIMS(self)[nd];
-            }
-            if (index[nd] < 0 ||
-                index[nd] >= PyArray_DIMS(self)[nd]) {
-                PyErr_SetString(PyExc_ValueError,
-                                "index out of bounds");
+
+        /* Convert the flat index into a multi-index */
+        for (idim = ndim-1; idim >= 0; --idim) {
+            multi_index[idim] = value % shape[idim];
+            value /= shape[idim];
+        }
+    }
+    /* A multi-index tuple */
+    else if (n == ndim) {
+        npy_intp value;
+
+        for (idim = 0; idim < ndim; ++idim) {
+            value = PyArray_PyIntAsIntp(PyTuple_GET_ITEM(args, idim));
+            if (value == -1 && PyErr_Occurred()) {
                 return NULL;
             }
-            loc += PyArray_STRIDES(self)[nd]*index[nd];
+            multi_index[idim] = value;
         }
-        ret = PyArray_DESCR(self)->f->setitem(obj, PyArray_DATA(self) + loc, self);
+    }
+    else {
+        PyErr_SetString(PyExc_ValueError,
+                        "incorrect number of indices for array");
+        return NULL;
     }
 
- finish:
-    if (ret < 0) {
+    if (PyArray_MultiIndexSetItem(self, multi_index, obj) < 0) {
         return NULL;
     }
-    Py_INCREF(Py_None);
-    return Py_None;
+    else {
+        Py_INCREF(Py_None);
+        return Py_None;
+    }
 }
 
 /* Sets the array values from another array as if they were flat */
@@ -1476,7 +1591,7 @@ array_setstate(PyArrayObject *self, PyObject *args)
     PyObject *rawdata = NULL;
     char *datastr;
     Py_ssize_t len;
-    intp size, dimensions[MAX_DIMS];
+    npy_intp size, dimensions[MAX_DIMS];
     int nd;
 
     PyArrayObject_fieldaccess *fa = (PyArrayObject_fieldaccess *)self;
@@ -1592,7 +1707,7 @@ array_setstate(PyArrayObject *self, PyObject *args)
         fa->dimensions = PyDimMem_NEW(3*nd);
         fa->strides = PyArray_DIMS(self) + nd;
         fa->maskna_strides = PyArray_DIMS(self) + 2*nd;
-        memcpy(PyArray_DIMS(self), dimensions, sizeof(intp)*nd);
+        memcpy(PyArray_DIMS(self), dimensions, sizeof(npy_intp)*nd);
         _array_fill_strides(PyArray_STRIDES(self), dimensions, nd,
                                PyArray_DESCR(self)->elsize,
                                (is_f_order ? NPY_ARRAY_F_CONTIGUOUS :
@@ -1610,7 +1725,7 @@ array_setstate(PyArrayObject *self, PyObject *args)
         /* Bytes are never interned */
         if (!_IsAligned(self) || swap) {
 #endif
-            intp num = PyArray_NBYTES(self);
+            npy_intp num = PyArray_NBYTES(self);
             fa->data = PyDataMem_NEW(num);
             if (PyArray_DATA(self) == NULL) {
                 fa->nd = 0;
@@ -1619,7 +1734,7 @@ array_setstate(PyArrayObject *self, PyObject *args)
                 return PyErr_NoMemory();
             }
             if (swap) { /* byte-swap on pickle-read */
-                intp numels = num / PyArray_DESCR(self)->elsize;
+                npy_intp numels = num / PyArray_DESCR(self)->elsize;
                 PyArray_DESCR(self)->f->copyswapn(PyArray_DATA(self),
                                         PyArray_DESCR(self)->elsize,
                                         datastr, PyArray_DESCR(self)->elsize,
diff --git a/numpy/core/src/multiarray/nditer_pywrap.c b/numpy/core/src/multiarray/nditer_pywrap.c
index 37aa8338b..b2ce78c11 100644
--- a/numpy/core/src/multiarray/nditer_pywrap.c
+++ b/numpy/core/src/multiarray/nditer_pywrap.c
@@ -2045,6 +2045,7 @@ npyiter_seq_item(NewNpyArrayIterObject *self, Py_ssize_t i)
     char *dataptr;
     PyArray_Descr *dtype;
     int has_external_loop;
+    Py_ssize_t i_orig = i;
 
     if (self->iter == NULL || self->finished) {
         PyErr_SetString(PyExc_ValueError,
@@ -2064,9 +2065,15 @@ npyiter_seq_item(NewNpyArrayIterObject *self, Py_ssize_t i)
      * before the first MASKNA operand.
      */
     nop = NpyIter_GetFirstMaskNAOp(self->iter);
+
+    /* Negative indexing */
+    if (i < 0) {
+        i += nop;
+    }
+
     if (i < 0 || i >= nop) {
         PyErr_Format(PyExc_IndexError,
-                "Iterator operand index %d is out of bounds", (int)i);
+                "Iterator operand index %d is out of bounds", (int)i_orig);
         return NULL;
     }
 
@@ -2113,8 +2120,6 @@ npyiter_seq_item(NewNpyArrayIterObject *self, Py_ssize_t i)
         return NULL;
     }
 
-    PyArray_UpdateFlags(ret, NPY_ARRAY_UPDATE_ALL);
-
     /* If this is a USE_MASKNA operand, include the mask */
     if (maskna_indices[i] >= 0) {
         PyArrayObject_fieldaccess *fret = (PyArrayObject_fieldaccess *)ret;
@@ -2131,6 +2136,8 @@ npyiter_seq_item(NewNpyArrayIterObject *self, Py_ssize_t i)
         fret->flags &= ~NPY_ARRAY_OWNMASKNA;
     }
 
+    PyArray_UpdateFlags(ret, NPY_ARRAY_UPDATE_ALL);
+
     return (PyObject *)ret;
 }
 
@@ -2198,6 +2205,8 @@ npyiter_seq_ass_item(NewNpyArrayIterObject *self, Py_ssize_t i, PyObject *v)
     PyArray_Descr *dtype;
     PyArrayObject *tmp;
     int ret, has_external_loop;
+    Py_ssize_t i_orig = i;
+
 
     if (v == NULL) {
         PyErr_SetString(PyExc_ValueError,
@@ -2223,14 +2232,20 @@ npyiter_seq_ass_item(NewNpyArrayIterObject *self, Py_ssize_t i, PyObject *v)
      * before the first MASKNA operand.
      */
     nop = NpyIter_GetFirstMaskNAOp(self->iter);
+
+    /* Negative indexing */
+    if (i < 0) {
+        i += nop;
+    }
+
     if (i < 0 || i >= nop) {
         PyErr_Format(PyExc_IndexError,
-                "Iterator operand index %d is out of bounds", (int)i);
+                "Iterator operand index %d is out of bounds", (int)i_orig);
         return -1;
     }
     if (!self->writeflags[i]) {
         PyErr_Format(PyExc_RuntimeError,
-                "Iterator operand %d is not writeable", (int)i);
+                "Iterator operand %d is not writeable", (int)i_orig);
         return -1;
     }
 
@@ -2276,7 +2291,9 @@ npyiter_seq_ass_item(NewNpyArrayIterObject *self, Py_ssize_t i, PyObject *v)
         ftmp->flags |= NPY_ARRAY_MASKNA;
         ftmp->flags &= ~NPY_ARRAY_OWNMASKNA;
     }
+
     PyArray_UpdateFlags(tmp, NPY_ARRAY_UPDATE_ALL);
+
     ret = PyArray_CopyObject(tmp, v);
     Py_DECREF(tmp);
     return ret;
diff --git a/numpy/core/tests/test_maskna.py b/numpy/core/tests/test_maskna.py
index 4a04123e7..3c188de21 100644
--- a/numpy/core/tests/test_maskna.py
+++ b/numpy/core/tests/test_maskna.py
@@ -73,6 +73,15 @@ def test_array_maskna_construction():
     assert_(a.flags.maskna)
     assert_equal(np.isna(a), True)
 
+def test_array_maskna_repr():
+    # Test some simple reprs with NA in them
+    a = np.array(np.NA, maskna=True)
+    assert_equal(repr(a), 'array(NA, maskna=True, dtype=float64)')
+    a = np.array([np.NA, 3], maskna=True)
+    assert_equal(repr(a), 'array([NA, 3], maskna=True)')
+    a = np.array([3.5, np.NA], maskna=True)
+    assert_equal(repr(a), 'array([ 3.5, NA], maskna=True)')
+
 def test_isna():
     # Objects which are not np.NA or ndarray all return False
     assert_equal(np.isna(True), False)
@@ -86,6 +95,50 @@ def test_isna():
     assert_equal(np.isna(np.NA(dtype='f4')), True)
     assert_equal(np.isna(np.NA(12,dtype='f4')), True)
 
+def test_array_maskna_item():
+    # With a zero-dimensional array
+    a = np.array(np.NA, maskna=True)
+
+    # Should return NA as the item
+    assert_equal(type(a.item()), np.NAType)
+
+    # Should be able to set the item
+    a.itemset(1.5)
+    assert_(not np.isna(a))
+    assert_equal(a, 1.5)
+    a.itemset(np.NA)
+    assert_(np.isna(a))
+
+    # With a one-dimensional array
+    a = np.array([1, np.NA, 2, np.NA], maskna=True)
+
+    # Should return the scalar or NA as the item
+    assert_(not np.isna(a.item(0)))
+    assert_equal(type(a.item(1)), np.NAType)
+
+    # Should be able to set the items
+    a.itemset(0, np.NA)
+    assert_(np.isna(a[0]))
+    a.itemset(1, 12)
+    assert_(not np.isna(a[1]))
+    assert_equal(a[1], 12)
+
+    # With a two-dimensional array
+    a = np.arange(6, maskna=True).reshape(2,3)
+    a[0,1] = np.NA
+    # Should return the scalar or NA as the item
+    assert_(not np.isna(a.item((0,0))))
+    assert_equal(type(a.item((0,1))), np.NAType)
+
+    # Should be able to set the items
+    a.itemset((0,1), 8)
+    assert_(not np.isna(a[0,1]))
+    assert_equal(a[0,1], 8)
+    a.itemset((1,1), np.NA)
+    assert_(np.isna(a[1,1]))
+
+
+
 def test_array_maskna_payload():
     # Single numbered index
     a = np.zeros((2,), maskna=True)
@@ -343,6 +396,10 @@ def test_array_maskna_reshape():
     assert_(not b.flags.ownmaskna)
     assert_equal(np.isna(b), [[0,0,0],[1,0,1]])
 
+    # Add a new axis using 'newaxis'
+    a = np.array(np.NA, maskna=True)
+    assert_equal(np.isna(a[np.newaxis]), [True])
+
 def test_array_maskna_view_NA_assignment_1D():
     a = np.arange(10)
     a_ref = a.copy()
@@ -788,6 +845,7 @@ def test_array_maskna_sum_prod_methods():
     assert_equal(res[~np.isna(res)], [4*8*7])
 
 def test_array_maskna_std_mean_methods():
+    return
     # ndarray.std, ndarray.mean
     a = np.array([[2, np.NA, 10],
                   [4, 8, 7],