#!/usr/bin/env pytest
###############################################################################
# $Id$
#
# Project:  GDAL/OGR Test Suite
# Purpose:  Test multidimensional support in MEM driver
# Author:   Even Rouault <even.rouault@spatialys.com>
#
###############################################################################
# Copyright (c) 2019, Even Rouault <even.rouault@spatialys.com>
#
# Permission is hereby granted, free of charge, to any person obtaining a
# copy of this software and associated documentation files (the "Software"),
# to deal in the Software without restriction, including without limitation
# the rights to use, copy, modify, merge, publish, distribute, sublicense,
# and/or sell copies of the Software, and to permit persons to whom the
# Software is furnished to do so, subject to the following conditions:
#
# The above copyright notice and this permission notice shall be included
# in all copies or substantial portions of the Software.
#
# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
# OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
# THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
# FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
# DEALINGS IN THE SOFTWARE.
###############################################################################

import array
import math
import struct
import sys

import gdaltest
import pytest

from osgeo import gdal, osr


###############################################################################
@pytest.fixture(autouse=True, scope="module")
def module_disable_exceptions():
    with gdaltest.disable_exceptions():
        yield


def test_mem_md_basic():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    assert ds
    assert ds.GetDescription() == "myds"

    rg = ds.GetRootGroup()
    assert rg
    assert rg.GetName() == "/"
    assert rg.GetFullName() == "/"
    assert not rg.GetMDArrayNames()
    assert not rg.GetGroupNames()
    assert not rg.GetAttributes()
    assert not rg.GetDimensions()
    assert not rg.OpenMDArray("not existing")
    assert not rg.OpenMDArrayFromFullname("not existing")
    assert not rg.OpenGroup("not existing")
    assert not rg.OpenGroupFromFullname("not existing")
    assert not rg.GetAttribute("not existing")
    assert not rg.GetVectorLayerNames()
    assert not rg.OpenVectorLayer("not existing")
    with gdal.ExceptionMgr(), pytest.raises(Exception):
        rg.OpenMDArray("not existing")
    with gdal.ExceptionMgr(), pytest.raises(Exception):
        rg.OpenMDArrayFromFullname("not existing")
    with gdal.ExceptionMgr(), pytest.raises(Exception):
        rg.OpenGroup("not existing")
    with gdal.ExceptionMgr(), pytest.raises(Exception):
        rg.OpenGroupFromFullname("not existing")
    with gdal.ExceptionMgr(), pytest.raises(Exception):
        rg.GetAttribute("not existing")
    with gdal.ExceptionMgr(), pytest.raises(Exception):
        rg.OpenVectorLayer("not existing")


def test_mem_md_subgroup():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    with gdal.quiet_errors():
        assert not rg.CreateGroup("")  # unnamed group not supported
    with pytest.raises(ValueError):
        assert not rg.CreateGroup(None)

    subg = rg.CreateGroup("subgroup")
    assert subg
    assert subg.GetName() == "subgroup"
    assert subg.GetFullName() == "/subgroup"
    assert rg.GetGroupNames() == ["subgroup"]
    assert rg.OpenGroup("subgroup").GetName() == "subgroup"

    subsubg = subg.CreateGroup("subsubgroup")
    assert subsubg.GetFullName() == "/subgroup/subsubgroup"

    subsubg = rg.OpenGroupFromFullname("/subgroup/subsubgroup")
    assert subsubg is not None
    assert subsubg.GetFullName() == "/subgroup/subsubgroup"

    subg.CreateMDArray("myarray", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte))
    array = rg.OpenMDArrayFromFullname("/subgroup/myarray")
    assert array is not None
    assert array.GetFullName() == "/subgroup/myarray"
    with gdal.ExceptionMgr(), pytest.raises(Exception):
        array.GetAttribute("not existing")

    copy_ds = drv.CreateCopy("", ds)
    assert copy_ds
    copy_rg = copy_ds.GetRootGroup()
    assert copy_rg
    assert copy_rg.GetGroupNames() == ["subgroup"]


def test_mem_md_array_unnamed_array():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    edt = gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    with gdal.quiet_errors():
        assert not rg.CreateMDArray("", [], edt)


def test_mem_md_array_duplicated_array_name():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    assert rg.CreateMDArray(
        "same_name", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    with gdal.quiet_errors():
        assert not rg.CreateMDArray(
            "same_name", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
        )


def test_mem_md_array_nodim():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    myarray = rg.CreateMDArray(
        "myarray", [], gdal.ExtendedDataType.Create(gdal.GDT_UInt16)
    )
    assert myarray
    assert myarray.GetName() == "myarray"
    assert myarray.GetFullName() == "/myarray"
    assert rg.GetMDArrayNames() == ["myarray"]
    assert rg.OpenMDArray("myarray")
    assert myarray.GetDimensionCount() == 0
    assert myarray.GetTotalElementsCount() == 1
    assert not myarray.GetDimensions()
    assert myarray.shape is None
    assert myarray.GetDataType().GetClass() == gdal.GEDTC_NUMERIC
    assert myarray.GetDataType().GetNumericDataType() == gdal.GDT_UInt16
    got_data = myarray.Read()
    assert len(got_data) == 2
    assert struct.unpack("H", got_data) == (0,)
    assert myarray.Write(struct.pack("H", 65535)) == gdal.CE_None
    got_data = myarray.Read()
    assert len(got_data) == 2
    assert struct.unpack("H", got_data) == (65535,)

    assert myarray.AdviseRead() == gdal.CE_None

    copy_ds = drv.CreateCopy("", ds)
    assert copy_ds
    copy_rg = copy_ds.GetRootGroup()
    assert copy_rg
    copy_myarray = copy_rg.OpenMDArray("myarray")
    assert copy_myarray
    assert copy_myarray.Read() == got_data

    assert len(myarray.GetCoordinateVariables()) == 0


def test_mem_md_array_single_dim():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim = rg.CreateDimension("dim0", "unspecified type", "unspecified direction", 2)
    dims = rg.GetDimensions()
    assert len(dims) == 1
    myarray = rg.CreateMDArray(
        "myarray", [dim], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    assert myarray
    assert myarray.GetName() == "myarray"
    assert rg.GetMDArrayNames() == ["myarray"]
    assert rg.OpenMDArray("myarray")
    assert myarray.GetDimensionCount() == 1
    assert myarray.GetTotalElementsCount() == 2
    got_dims = myarray.GetDimensions()
    assert len(got_dims) == 1
    assert myarray.shape == (2,)
    assert got_dims[0].GetName() == "dim0"
    assert got_dims[0].GetType() == "unspecified type"
    assert got_dims[0].GetDirection() == "unspecified direction"
    assert got_dims[0].GetSize() == 2
    got_data = myarray.Read()
    assert len(got_data) == 2
    assert struct.unpack("B" * 2, got_data) == (0, 0)

    # Test writing a list with only integer values
    assert myarray.Write([0, 1]) == gdal.CE_None
    assert struct.unpack("B" * 2, myarray.Read()) == (0, 1)

    # Unsupported type in list
    with pytest.raises(Exception):
        myarray.Write([0, "aa"])

    # Test writing a list with a mix of int and double
    assert myarray.Write([0, 2.1]) == gdal.CE_None
    assert struct.unpack("B" * 2, myarray.Read()) == (0, 2)

    # Test writing a list with only double values
    assert myarray.Write([1.0, 2.1]) == gdal.CE_None
    assert struct.unpack("B" * 2, myarray.Read()) == (1, 2)

    # Test writing a list with only integer values, but one out of int32 range
    assert myarray.Write([1, 1 << 31]) == gdal.CE_None
    assert struct.unpack("B" * 2, myarray.Read()) == (1, 255)

    # Test writing a array
    for typecode, in_ar, out_tuple in [
        ("B", [1, 2], (1, 2)),
        ("b", [-128, 127], (0, 127)),
        ("h", [-32768, 32767], (0, 255)),
        ("H", [0, 65535], (0, 255)),
        ("i", [-(1 << 31), (1 << 31) - 1], (0, 255)),
        ("I", [0, (1 << 32) - 1], (0, 255)),
        ("f", [1.0, 2.1], (1, 2)),
        ("d", [2.0, 3.1], (2, 3)),
    ]:
        assert myarray.Write(array.array(typecode, in_ar)) == gdal.CE_None
        assert struct.unpack("B" * 2, myarray.Read()) == out_tuple

    assert myarray.AdviseRead() == gdal.CE_None

    attr = myarray.CreateAttribute(
        "attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    assert attr
    assert attr.GetFullName() == "/myarray/attr"

    assert myarray.GetNoDataValueAsDouble() is None
    assert myarray.SetNoDataValueDouble(1) == gdal.CE_None
    assert myarray.GetNoDataValueAsDouble() == 1
    assert myarray.GetNoDataValueAsInt64() == 1
    assert myarray.GetNoDataValueAsUInt64() == 1
    assert myarray.SetNoDataValue(1) == gdal.CE_None
    assert myarray.GetNoDataValue() == 1
    assert myarray.SetNoDataValueRaw(struct.pack("B", 127)) == gdal.CE_None
    with gdal.quiet_errors():
        assert myarray.SetNoDataValueRaw(struct.pack("h", 127)) != gdal.CE_None
    assert struct.unpack("B", myarray.GetNoDataValueAsRaw()) == (127,)

    assert myarray.GetScale() is None
    assert myarray.GetOffset() is None
    assert myarray.GetScaleStorageType() == gdal.GDT_Unknown
    assert myarray.GetOffsetStorageType() == gdal.GDT_Unknown

    assert myarray.SetScale(2.5) == gdal.CE_None
    assert myarray.GetScale() == 2.5
    assert myarray.GetScaleStorageType() == gdal.GDT_Unknown
    assert myarray.SetScale(2.5, storageType=gdal.GDT_Float32) == gdal.CE_None
    assert myarray.GetScaleStorageType() == gdal.GDT_Float32

    assert myarray.SetOffset(1.5) == gdal.CE_None
    assert myarray.GetOffset() == 1.5
    assert myarray.GetOffsetStorageType() == gdal.GDT_Unknown
    assert myarray.SetOffset(1.5, storageType=gdal.GDT_Float32) == gdal.CE_None
    assert myarray.GetOffsetStorageType() == gdal.GDT_Float32

    def my_cbk(pct, _, arg):
        assert pct >= tab[0]
        tab[0] = pct
        return 1

    got_data = myarray.Read()

    tab = [0]
    copy_ds = drv.CreateCopy("", ds, callback=my_cbk, callback_data=tab)
    assert tab[0] == 1
    assert copy_ds
    copy_rg = copy_ds.GetRootGroup()
    assert copy_rg
    copy_myarray = copy_rg.OpenMDArray("myarray")
    assert copy_myarray
    assert copy_myarray.Read() == got_data
    assert copy_myarray.GetNoDataValueAsRaw() == myarray.GetNoDataValueAsRaw()

    assert myarray.DeleteNoDataValue() == gdal.CE_None
    assert myarray.GetNoDataValueAsDouble() is None
    with gdal.quiet_errors():
        assert myarray.GetNoDataValueAsString() is None

    assert myarray.SetUnit("foo") == gdal.CE_None
    assert myarray.GetUnit() == "foo"


def test_mem_md_array_string():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim = rg.CreateDimension("dim0", "unspecified type", "unspecified direction", 2)
    var = rg.CreateMDArray("var", [dim], gdal.ExtendedDataType.CreateString())
    assert var
    assert var.Read() == [None, None]
    assert var.Write(["", "0123456789"]) == gdal.CE_None
    assert var.GetNoDataValueAsString() is None
    assert var.SetNoDataValueString(None) == gdal.CE_None
    assert var.GetNoDataValueAsString() is None
    assert var.SetNoDataValueString("123") == gdal.CE_None
    assert var.GetNoDataValueAsString() == "123"
    var = rg.OpenMDArray("var")
    assert var
    assert var.Read() == ["", "0123456789"]


def test_mem_md_datatypes():

    dt_byte = gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    assert dt_byte.GetClass() == gdal.GEDTC_NUMERIC
    assert dt_byte.GetName() == ""
    assert dt_byte.GetNumericDataType() == gdal.GDT_Byte
    assert dt_byte.GetSize() == 1
    assert dt_byte.CanConvertTo(dt_byte)
    with pytest.raises(ValueError):
        assert dt_byte.CanConvertTo(None)
    assert dt_byte == gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    assert not dt_byte != gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    assert dt_byte.Equals(dt_byte)
    with pytest.raises(ValueError):
        assert dt_byte.Equals(None)
    assert not dt_byte.GetComponents()

    dt_cint32 = gdal.ExtendedDataType.Create(gdal.GDT_CInt32)
    assert dt_cint32.GetClass() == gdal.GEDTC_NUMERIC
    assert dt_cint32.GetName() == ""
    assert dt_cint32.GetNumericDataType() == gdal.GDT_CInt32
    assert dt_cint32.GetSize() == 2 * 4
    assert dt_cint32.CanConvertTo(dt_cint32)
    assert dt_cint32.Equals(dt_cint32)
    assert dt_cint32.CanConvertTo(dt_byte)
    assert not dt_cint32.Equals(dt_byte)
    assert not dt_cint32.GetComponents()

    dt_string = gdal.ExtendedDataType.CreateString()
    assert dt_string.GetClass() == gdal.GEDTC_STRING
    assert dt_string.GetName() == ""
    assert dt_string.GetNumericDataType() == gdal.GDT_Unknown
    assert dt_string.GetSize() in (4, 8)  # depends on 32 vs 64 build
    assert dt_string.GetMaxStringLength() == 0
    assert dt_string.CanConvertTo(dt_string)
    assert dt_string.Equals(dt_string)
    assert not dt_string.Equals(dt_byte)
    assert not dt_byte.Equals(dt_string)
    assert not dt_string.GetComponents()

    dt_string_limited_size = gdal.ExtendedDataType.CreateString(10)
    assert dt_string_limited_size.GetMaxStringLength() == 10

    comp0 = gdal.EDTComponent.Create(
        "x", 0, gdal.ExtendedDataType.Create(gdal.GDT_Int16)
    )
    comp1 = gdal.EDTComponent.Create(
        "y", 4, gdal.ExtendedDataType.Create(gdal.GDT_Int32)
    )

    with gdal.quiet_errors():
        assert gdal.ExtendedDataType.CreateCompound("mytype", 8, []) is None
        assert (
            gdal.ExtendedDataType.CreateCompound("mytype", 2000 * 1000 * 1000, [comp0])
            is None
        )

    compound_dt = gdal.ExtendedDataType.CreateCompound("mytype", 8, [comp0, comp1])
    assert compound_dt.GetClass() == gdal.GEDTC_COMPOUND
    assert compound_dt.GetName() == "mytype"
    assert compound_dt.GetNumericDataType() == gdal.GDT_Unknown
    assert compound_dt.GetSize() == 8
    comps = compound_dt.GetComponents()
    assert len(comps) == 2
    assert comps[0].GetName() == "x"
    assert comps[0].GetOffset() == 0
    assert comps[0].GetType().GetNumericDataType() == gdal.GDT_Int16
    assert comps[1].GetName() == "y"
    assert comps[1].GetOffset() == 4
    assert comps[1].GetType().GetNumericDataType() == gdal.GDT_Int32
    assert compound_dt.CanConvertTo(compound_dt)
    assert compound_dt.Equals(compound_dt)
    assert not compound_dt.Equals(dt_byte)
    assert not dt_byte.Equals(compound_dt)

    with gdal.quiet_errors():
        # Too short size
        assert not gdal.ExtendedDataType.CreateCompound("mytype", 7, [comp0, comp1])

        # Too big size
        assert not gdal.ExtendedDataType.CreateCompound(
            "mytype", 1 << 30, [comp0, comp1]
        )

        # Wrongly ordered
        assert not gdal.ExtendedDataType.CreateCompound("mytype", 8, [comp1, comp0])

        # Empty
        assert not gdal.ExtendedDataType.CreateCompound("mytype", 0, [])

    other_compound_dt = gdal.ExtendedDataType.CreateCompound("mytype", 8, [comp0])
    assert not compound_dt.Equals(other_compound_dt)
    assert not other_compound_dt.Equals(compound_dt)
    assert compound_dt.CanConvertTo(other_compound_dt)
    assert not other_compound_dt.CanConvertTo(compound_dt)


def test_mem_md_array_compoundtype():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim = rg.CreateDimension("dim0", None, None, 2)
    comp0 = gdal.EDTComponent.Create(
        "x", 0, gdal.ExtendedDataType.Create(gdal.GDT_Int16)
    )
    comp1 = gdal.EDTComponent.Create(
        "y", 4, gdal.ExtendedDataType.Create(gdal.GDT_Int32)
    )
    dt = gdal.ExtendedDataType.CreateCompound("mytype", 8, [comp0, comp1])
    myarray = rg.CreateMDArray("myarray", [dim], dt)
    assert myarray
    dt = myarray.GetDataType()
    assert dt.GetName() == "mytype"
    assert dt.GetClass() == gdal.GEDTC_COMPOUND
    assert dt.GetSize() == 8
    comps = dt.GetComponents()
    assert len(comps) == 2
    assert comps[0].GetName() == "x"
    assert comps[0].GetOffset() == 0
    assert comps[0].GetType().GetNumericDataType() == gdal.GDT_Int16
    assert comps[1].GetName() == "y"
    assert comps[1].GetOffset() == 4
    assert comps[1].GetType().GetNumericDataType() == gdal.GDT_Int32

    assert (
        myarray.Write(struct.pack("hi" * 2, 32767, 1000000, -32768, -1000000))
        == gdal.CE_None
    )
    assert struct.unpack("hi" * 2, myarray.Read()) == (32767, 1000000, -32768, -1000000)

    extract_compound_dt = gdal.ExtendedDataType.CreateCompound(
        "mytype",
        4,
        [
            gdal.EDTComponent.Create(
                "y", 0, gdal.ExtendedDataType.Create(gdal.GDT_Int32)
            )
        ],
    )
    got_data = myarray.Read(buffer_datatype=extract_compound_dt)
    assert len(got_data) == 8
    assert struct.unpack("i" * 2, got_data) == (1000000, -1000000)

    with gdal.quiet_errors():
        assert not myarray.GetView('["z')
        assert not myarray.GetView('["z"')
        assert not myarray.GetView('["z"]')
    y_ar = myarray["y"]
    assert y_ar
    assert y_ar.GetDimensionCount() == 1
    assert y_ar.GetDataType().GetClass() == gdal.GEDTC_NUMERIC
    assert y_ar.GetDataType().GetNumericDataType() == gdal.GDT_Int32
    assert y_ar.GetBlockSize() == myarray.GetBlockSize()

    assert y_ar.GetSpatialRef() == myarray.GetSpatialRef()
    assert y_ar.GetNoDataValueAsRaw() == myarray.GetNoDataValueAsRaw()
    assert y_ar.GetUnit() == myarray.GetUnit()
    assert y_ar.GetScale() == myarray.GetScale()
    assert y_ar.GetOffset() == myarray.GetOffset()

    assert myarray.SetUnit("foo") == gdal.CE_None
    assert myarray.SetScale(1) == gdal.CE_None
    assert myarray.SetOffset(2) == gdal.CE_None
    assert myarray.SetNoDataValueRaw(struct.pack("hi", 32767, 1000000)) == gdal.CE_None
    sr = osr.SpatialReference()
    sr.ImportFromEPSG(32631)
    assert myarray.SetSpatialRef(sr) == gdal.CE_None

    assert y_ar.GetUnit() == myarray.GetUnit()
    assert y_ar.GetScale() == myarray.GetScale()
    assert y_ar.GetOffset() == myarray.GetOffset()
    assert y_ar.GetNoDataValueAsRaw() == struct.pack("i", 1000000)
    assert y_ar.GetSpatialRef().IsSame(myarray.GetSpatialRef())

    got_data = y_ar.Read()
    assert len(got_data) == 8
    assert struct.unpack("i" * 2, got_data) == (1000000, -1000000)

    with gdal.quiet_errors():
        assert not y_ar.GetView('["y"]')

    assert y_ar.AdviseRead() == gdal.CE_None

    y_ar = myarray["y"][1]
    got_data = y_ar.Read()
    assert len(got_data) == 4
    assert struct.unpack("i" * 1, got_data) == (-1000000,)

    y_ar = myarray[1]["y"]
    got_data = y_ar.Read()
    assert len(got_data) == 4
    assert struct.unpack("i" * 1, got_data) == (-1000000,)


def test_mem_md_array_3_dim():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim0 = rg.CreateDimension("dim0", None, None, 2)
    dim1 = rg.CreateDimension("dim1", None, None, 3)
    dim2 = rg.CreateDimension("dim2", None, None, 4)
    myarray = rg.CreateMDArray(
        "myarray", [dim0, dim1, dim2], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    assert myarray
    assert myarray.GetName() == "myarray"
    assert rg.GetMDArrayNames() == ["myarray"]
    assert rg.OpenMDArray("myarray")
    assert myarray.GetDimensionCount() == 3
    assert myarray.GetTotalElementsCount() == 24
    assert myarray.shape == (2, 3, 4)

    assert not myarray.GetSpatialRef()
    assert myarray.SetSpatialRef(None) == gdal.CE_None
    sr = osr.SpatialReference()
    sr.ImportFromEPSG(32631)
    assert myarray.SetSpatialRef(sr) == gdal.CE_None
    assert myarray.GetSpatialRef() is not None

    data = array.array("B", list(range(24))).tobytes()
    assert myarray.Write(data) == gdal.CE_None

    got_data = myarray.Read()
    assert got_data == data

    got_data = myarray.Read(
        [0, 0, 0],
        [2, 3, 4],
        [1, 1, 1],
        [12, 4, 1],
        gdal.ExtendedDataType.Create(gdal.GDT_Byte),
    )
    assert got_data == data

    data = myarray.Read(
        array_start_idx=[1, 0, 2], count=[1, 2, 2], array_step=[1, 2, 1]
    )
    assert struct.unpack("B" * len(data), data) == (14, 15, 22, 23)

    data = myarray.Read(
        array_start_idx=[1, 0, 3], count=[1, 2, 2], array_step=[1, 2, -1]
    )
    assert struct.unpack("B" * len(data), data) == (15, 14, 23, 22)

    data = struct.pack("d" * 1, 25.0)
    float64dt = gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    assert (
        myarray.Write(data, [1, 2, 3], [1, 1, 1], buffer_datatype=float64dt)
        == gdal.CE_None
    )

    got_data = myarray.Read([1, 2, 3], [1, 1, 1], buffer_datatype=float64dt)
    got_data = struct.unpack("d" * 1, got_data)[0]
    assert got_data == 25.0


def test_mem_md_array_4_dim():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim0 = rg.CreateDimension("dim0", None, None, 2)
    dim1 = rg.CreateDimension("dim1", None, None, 3)
    dim2 = rg.CreateDimension("dim2", None, None, 4)
    dim3 = rg.CreateDimension("dim3", None, None, 5)
    myarray = rg.CreateMDArray(
        "myarray",
        [dim0, dim1, dim2, dim3],
        gdal.ExtendedDataType.Create(gdal.GDT_Int16),
    )
    assert myarray
    assert myarray.GetName() == "myarray"
    assert rg.GetMDArrayNames() == ["myarray"]
    assert rg.OpenMDArray("myarray")
    assert myarray.GetDimensionCount() == 4
    assert myarray.GetTotalElementsCount() == 2 * 3 * 4 * 5

    data = array.array("h", [-i for i in range(2 * 3 * 4 * 5)]).tobytes()
    assert myarray.Write(data) == gdal.CE_None

    got_data = myarray.Read()
    assert got_data == data

    data = myarray.Read(
        array_start_idx=[1, 0, 2, 3], count=[1, 2, 2, 2], array_step=[1, 2, 1, 1]
    )
    assert struct.unpack("h" * (len(data) // 2), data) == (
        -73,
        -74,
        -78,
        -79,
        -113,
        -114,
        -118,
        -119,
    )


def test_mem_md_copy_array():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim0 = rg.CreateDimension("dim0", None, None, 20)
    dim1 = rg.CreateDimension("dim1", None, None, 30)
    dim2 = rg.CreateDimension("dim2", None, None, 10)
    dim3 = rg.CreateDimension("dim3", None, None, 51)
    myarray = rg.CreateMDArray(
        "myarray",
        [dim0, dim1, dim2, dim3],
        gdal.ExtendedDataType.Create(gdal.GDT_UInt32),
    )

    data = array.array("I", list(range(myarray.GetTotalElementsCount()))).tobytes()
    assert myarray.Write(data) == gdal.CE_None

    def my_cbk(pct, _, arg):
        assert pct > tab[0]
        tab[0] = pct
        return 1

    tab = [0]
    with gdaltest.config_option("GDAL_SWATH_SIZE", str(100 * 1000)):
        copy_ds = drv.CreateCopy("", ds, callback=my_cbk, callback_data=tab)
    assert tab[0] == 1
    assert copy_ds
    copy_rg = copy_ds.GetRootGroup()
    assert copy_rg
    copy_myarray = copy_rg.OpenMDArray("myarray")
    assert copy_myarray
    assert copy_myarray.Read() == data


def test_mem_md_array_read_write_errors():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim0 = rg.CreateDimension("dim0", None, None, 2)
    dim1 = rg.CreateDimension("dim1", None, None, 3)
    dim2 = rg.CreateDimension("dim2", None, None, 4)
    myarray = rg.CreateMDArray(
        "myarray", [dim0, dim1, dim2], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    assert myarray

    assert myarray.Read([0, 0, 0], [1, 1, 1], [1, 1, 1], [0, 0, 0])
    with gdal.quiet_errors():

        # Invalid number of values in array_idx array
        assert not myarray.Read([0, 0], [1, 1, 1], [1, 1, 1], [0, 0, 0])

        # Invalid number of values in count array
        assert not myarray.Read([0, 0, 0], [1, 1], [1, 1, 1], [0, 0, 0])

        # Invalid number of values in step array
        assert not myarray.Read([0, 0, 0], [1, 1, 1], [1, 1], [0, 0, 0])

        # Invalid number of values in buffer_stride array
        assert not myarray.Read([0, 0, 0], [1, 1, 1], [1, 1, 1], [0, 0])

        # Invalid start_index[0]
        assert not myarray.Read([-1, 0, 0], [1, 1, 1], [1, 1, 1], [0, 0, 0])

        # Invalid start_index[0]
        assert not myarray.Read([2, 0, 0], [1, 1, 1], [1, 1, 1], [0, 0, 0])

        # Invalid count[0]
        assert not myarray.Read([0, 0, 0], [-1, 1, 1], [1, 1, 1], [0, 0, 0])

        # Invalid count[0]
        assert not myarray.Read([0, 0, 0], [0, 1, 1], [1, 1, 1], [0, 0, 0])

        # Invalid count[0]
        assert not myarray.Read([0, 0, 0], [3, 1, 1], [1, 1, 1], [0, 0, 0])

        # start_idx[0] + (count[0]-1) * step[0] >= dim[0].size
        assert not myarray.Read([0, 0, 0], [2, 1, 1], [2, 1, 1], [0, 0, 0])
        assert not myarray.Read([0, 0, 0], [2 << 50, 1, 1], [2 << 50, 1, 1], [0, 0, 0])

        # Overflow with step[0]
        assert not myarray.Read([0, 0, 0], [2, 1, 1], [-(1 << 63), 1, 1], [0, 0, 0])

        # start_idx[0] + (count[0]-1) * step[0] < 0
        assert not myarray.Read([0, 0, 0], [2, 1, 1], [-1, 1, 1], [0, 0, 0])
        assert not myarray.Read([0, 0, 0], [2, 1, 1], [-2 << 50, 1, 1], [0, 0, 0])

        # Too big stride
        assert not myarray.Read([1, 0, 0], [2, 1, 1], [1, 1, 1], [(1 << 63) - 1, 0, 0])

        # Negative stride not supported in SWIG bindings
        assert not myarray.Read([0, 0, 0], [1, 1, 1], [1, 1, 1], [-1, 0, 0])

    data = struct.pack("d" * 1, 25.0)
    float64dt = gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    with gdal.quiet_errors():
        assert myarray.Write("", [1, 2, 3], [1, 1, 1]) == gdal.CE_Failure
        assert (
            myarray.Write(data[0:7], [1, 2, 3], [1, 1, 1], buffer_datatype=float64dt)
            == gdal.CE_Failure
        )


def test_mem_md_invalid_dims():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    assert rg.CreateDimension("dim1", None, None, 1)
    with gdal.quiet_errors():
        # empty name
        assert not rg.CreateDimension("", None, None, 1)
        # existing dim
        assert not rg.CreateDimension("dim1", None, None, 1)


def test_mem_md_array_invalid_args():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    edt = gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    dim = rg.CreateDimension("dim0", None, None, 1)
    with pytest.raises(TypeError):
        rg.CreateMDArray("myarray", None, edt)
    with pytest.raises((TypeError, SystemError)):
        rg.CreateMDArray("myarray", [None], edt)
    with pytest.raises((TypeError, SystemError)):
        rg.CreateMDArray("myarray", [1], edt)
    with pytest.raises(ValueError):
        rg.CreateMDArray("myarray", [dim], None)
    with pytest.raises(ValueError):
        rg.CreateMDArray(None, [dim], edt)


def test_mem_md_array_too_large():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim = rg.CreateDimension("dim0", None, None, (1 << 64) - 1)
    with gdal.quiet_errors():
        assert not rg.CreateMDArray(
            "myarray", [dim], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
        )


def test_mem_md_array_too_large_overflow_dim():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim0 = rg.CreateDimension("dim0", None, None, 1 << 25)
    dim1 = rg.CreateDimension("dim1", None, None, 1 << 25)
    dim2 = rg.CreateDimension("dim2", None, None, 1 << 25)
    with gdal.quiet_errors():
        assert not rg.CreateMDArray(
            "myarray", [dim0, dim1, dim2], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
        )


def test_mem_md_array_30dim():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dims = []
    for i in range(30):
        dim = rg.CreateDimension("dim%d" % i, None, None, 1 + (i % 2))
        dims.append(dim)
    myarray = rg.CreateMDArray(
        "myarray", dims, gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    assert myarray
    assert myarray.GetDimensionCount() == 30
    assert myarray.GetTotalElementsCount() == 2**15
    got_data = myarray.Read()
    assert len(got_data) == 2**15


def test_mem_md_array_32dim():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dims = []
    for i in range(32):
        dim = rg.CreateDimension("dim%d" % i, None, None, 1 + (i % 2))
        dims.append(dim)
    myarray = rg.CreateMDArray(
        "myarray", dims, gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    assert len(myarray.Read()) == 2**16


def test_mem_md_group_attribute_single_numeric():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    float64dt = gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    with gdal.quiet_errors():
        assert not rg.CreateAttribute("", [1], float64dt)  # unnamed attr not supported
    with pytest.raises(ValueError):
        rg.CreateAttribute(None, [1], float64dt)

    attr = rg.CreateAttribute("attr", [1], float64dt)
    assert attr
    assert attr.GetName() == "attr"
    assert rg.GetAttribute("attr").GetName() == "attr"
    assert rg.GetAttributes()[0].GetName() == "attr"
    assert attr.GetDimensionsSize() == [1]
    assert attr.GetTotalElementsCount() == 1
    assert attr.GetDataType().GetClass() == gdal.GEDTC_NUMERIC
    assert attr.GetDataType().GetNumericDataType() == gdal.GDT_Float64
    assert attr.Read() == 0.0
    assert attr.ReadAsString() == "0"
    assert attr.Write(1) == gdal.CE_None
    assert attr.Read() == 1
    assert attr.Write(1.25) == gdal.CE_None
    assert attr.Read() == 1.25
    assert attr.ReadAsDouble() == 1.25
    assert attr.ReadAsDoubleArray() == (1.25,)
    assert attr.ReadAsString() == "1.25"
    assert attr.ReadAsStringArray() == ["1.25"]
    assert attr.Write([2]) == gdal.CE_None
    assert attr.Read() == 2
    assert attr.Write([2.25]) == gdal.CE_None
    assert attr.Read() == 2.25
    with gdal.quiet_errors():
        assert attr.Write([]) != gdal.CE_None
        assert attr.Write([1, 2]) != gdal.CE_None


def test_mem_md_group_attribute_multiple_numeric():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    float64dt = gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    attr = rg.CreateAttribute("attr", [2, 3], float64dt)
    assert attr
    assert attr.GetFullName() == "/_GLOBAL_/attr"
    assert attr.GetDimensionsSize() == [2, 3]
    assert attr.GetTotalElementsCount() == 6
    assert attr.Read() == (0.0, 0.0, 0.0, 0.0, 0.0, 0.0)
    assert attr.Write([1.1, 2, 3, 4, 5, 6]) == gdal.CE_None
    assert attr.Read() == (1.1, 2, 3, 4, 5, 6)
    assert attr.ReadAsIntArray() == (1, 2, 3, 4, 5, 6)
    assert attr.ReadAsInt() == 1
    assert attr.ReadAsDouble() == 1.1

    subg = rg.CreateGroup("subgroup")
    attr = subg.CreateAttribute("attr", [2, 3], float64dt)
    assert attr
    assert attr.GetFullName() == "/subgroup/_GLOBAL_/attr"


def test_mem_md_group_attribute_single_string():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    attr = rg.CreateAttribute("attr", [1], gdal.ExtendedDataType.CreateString())
    assert attr
    assert attr.Read() is None
    assert attr.ReadAsStringArray() == [""]
    assert attr.Write("foo") == gdal.CE_None
    assert attr.Read() == "foo"
    assert attr.ReadAsStringArray() == ["foo"]
    assert attr.Write(["bar"]) == gdal.CE_None
    assert attr.Read() == "bar"


def test_mem_md_group_attribute_string_json():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    attr = rg.CreateAttribute(
        "attr", [], gdal.ExtendedDataType.CreateString(0, gdal.GEDTST_JSON)
    )
    assert attr
    assert attr.GetDataType().GetSubType() == gdal.GEDTST_JSON
    assert attr.Read() is None
    assert attr.Write({"foo": "bar"}) == gdal.CE_None
    assert attr.Read() == {"foo": "bar"}


def test_mem_md_group_attribute_multiple_string():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    attr = rg.CreateAttribute("attr", [2, 3], gdal.ExtendedDataType.CreateString())
    assert attr
    assert attr.GetDimensionsSize() == [2, 3]
    assert attr.GetTotalElementsCount() == 6
    assert attr.Read() == ["", "", "", "", "", ""]
    assert attr.Write(["foo", "bar", "baz", "FOO", "BAR", "BAZ"]) == gdal.CE_None
    assert attr.Read() == ["foo", "bar", "baz", "FOO", "BAR", "BAZ"]


def test_mem_md_array_attribute():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    myarray = rg.CreateMDArray(
        "myarray", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )

    float64dt = gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    with gdal.quiet_errors():
        assert not myarray.CreateAttribute(
            "", [1], float64dt
        )  # unnamed attr not supported
    with pytest.raises(ValueError):
        myarray.CreateAttribute(None, [1], float64dt)

    attr = myarray.CreateAttribute("attr", [1], float64dt)
    assert attr
    assert attr.GetName() == "attr"
    assert myarray.GetAttribute("attr").GetName() == "attr"
    assert myarray.GetAttributes()[0].GetName() == "attr"
    assert not myarray.GetAttribute("not existing")


def test_mem_md_array_slice():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim_2 = rg.CreateDimension("dim_2", None, None, 2)
    dim_3 = rg.CreateDimension("dim_3", None, None, 3)
    dim_4 = rg.CreateDimension("dim_4", None, None, 4)

    ar = rg.CreateMDArray("nodim", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte))
    assert ar.Write(struct.pack("B", 1)) == gdal.CE_None
    with gdal.quiet_errors():
        assert not ar[:]

    ar = rg.CreateMDArray(
        "array", [dim_2, dim_3, dim_4], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    data = array.array("B", list(range(24))).tobytes()
    assert ar.Write(data) == gdal.CE_None

    with pytest.raises(Exception):
        ar.GetView(None)

    attr = ar.CreateAttribute(
        "attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    assert attr.Write(1) == gdal.CE_None

    with gdal.quiet_errors():
        assert not ar.GetView("")
        assert not ar.GetView("x")
        assert not ar.GetView("[")
        assert not ar.GetView("[]")
        assert not ar.GetView("[foo]")
        assert not ar.GetView("[1,2,3,4]")
        assert not ar.GetView("[1,2,3,:]")
        assert not ar.GetView("[1:2:3:4]")
        assert not ar.GetView("[...,...]")
        assert not ar.GetView("[...,0,...]")
        assert not ar.GetView("[-3]")
        assert not ar.GetView("[2]")
        assert not ar.GetView("[::0]")
        assert not ar.GetView("[0:0:1]")
        assert not ar.GetView("[1:1:1]")
        assert not ar.GetView("[1:1:-1]")
        assert not ar.GetView("[0,0,0][:]")

    sliced_ar = ar[1]
    assert sliced_ar
    assert sliced_ar.GetAttribute("attr") is not None
    assert len(sliced_ar.GetAttributes()) == 1
    assert sliced_ar.GetBlockSize() == ar.GetBlockSize()[1:]
    assert sliced_ar.GetDimensionCount() == ar.GetDimensionCount() - 1
    assert sliced_ar.GetDataType() == ar.GetDataType()
    dims = sliced_ar.GetDimensions()
    assert dims[0].GetName() == "dim_3"
    assert dims[1].GetName() == "dim_4"
    assert sliced_ar.Read() == ar.Read(array_start_idx=[1, 0, 0], count=[1, 3, 4])

    sliced_ar = ar[1, 2, 3]
    assert sliced_ar.GetDimensionCount() == 0
    assert sliced_ar.Read() == ar.Read(array_start_idx=[1, 2, 3], count=[1, 1, 1])

    sliced_ar = ar[1, "newaxis", 2, 3]
    assert sliced_ar.GetDimensionCount() == 1
    assert sliced_ar.Read() == ar.Read(array_start_idx=[1, 2, 3], count=[1, 1, 1])

    try:
        import numpy as np

        sliced_ar = ar[1, np.newaxis, 2, 3]
        assert sliced_ar.GetDimensionCount() == 1
        assert sliced_ar.Read() == ar.Read(array_start_idx=[1, 2, 3], count=[1, 1, 1])
    except ImportError:
        pass

    sliced_ar = ar[1][2][3]
    assert sliced_ar.GetDimensionCount() == 0
    assert sliced_ar.Read() == ar.Read(array_start_idx=[1, 2, 3], count=[1, 1, 1])

    sliced_ar = ar.GetView("[1][2][3]")
    assert sliced_ar.GetDimensionCount() == 0
    assert sliced_ar.Read() == ar.Read(array_start_idx=[1, 2, 3], count=[1, 1, 1])

    orig_data = sliced_ar.Read()
    assert sliced_ar.Write(struct.pack("B", 123)) == gdal.CE_None
    assert struct.unpack("B", sliced_ar.Read()) == (123,)
    assert sliced_ar.Write(orig_data) == gdal.CE_None

    sliced_ar = ar[...]
    assert sliced_ar.GetDimensionCount() == ar.GetDimensionCount()
    assert sliced_ar.Read() == ar.Read()
    assert sliced_ar.Read(
        array_start_idx=[1, 0, 0], count=[1, 3, 2], array_step=[1, 1, 2]
    ) == ar.Read(array_start_idx=[1, 0, 0], count=[1, 3, 2], array_step=[1, 1, 2])

    assert sliced_ar.GetUnit() == ar.GetUnit()
    assert sliced_ar.GetScale() == ar.GetScale()
    assert sliced_ar.GetOffset() == ar.GetOffset()
    assert sliced_ar.GetNoDataValueAsRaw() == ar.GetNoDataValueAsRaw()
    assert sliced_ar.GetSpatialRef() is None

    assert ar.SetUnit("foo") == gdal.CE_None
    assert ar.SetScale(1) == gdal.CE_None
    assert ar.SetOffset(2) == gdal.CE_None
    assert ar.SetNoDataValueDouble(1) == gdal.CE_None
    sr = osr.SpatialReference()
    sr.ImportFromEPSG(32631)
    assert ar.SetSpatialRef(sr) == gdal.CE_None
    assert ar.GetSpatialRef()

    assert sliced_ar.GetUnit() == ar.GetUnit()
    assert sliced_ar.GetScale() == ar.GetScale()
    assert sliced_ar.GetOffset() == ar.GetOffset()
    assert sliced_ar.GetNoDataValueAsRaw() == ar.GetNoDataValueAsRaw()
    assert sliced_ar.GetSpatialRef().IsSame(ar.GetSpatialRef())

    sliced_ar = ar[1, ...]
    assert sliced_ar.Read() == ar.Read(array_start_idx=[1, 0, 0], count=[1, 3, 4])

    sliced_ar = ar[..., 1]
    assert sliced_ar.GetDimensionCount() == ar.GetDimensionCount() - 1
    assert sliced_ar.Read() == ar.Read(array_start_idx=[0, 0, 1], count=[2, 3, 1])

    sliced_ar = ar[1, ..., 2]
    assert sliced_ar.GetDimensionCount() == ar.GetDimensionCount() - 2
    assert sliced_ar.Read() == ar.Read(array_start_idx=[1, 0, 2], count=[1, 3, 1])

    sliced_ar = ar.GetView("[:]")
    assert sliced_ar.GetDimensionCount() == ar.GetDimensionCount()
    assert sliced_ar.Read() == ar.Read()

    sliced_ar = ar.GetView("[:,:,1]")
    assert sliced_ar.GetDimensionCount() == ar.GetDimensionCount() - 1
    assert sliced_ar.Read() == ar.Read(array_start_idx=[0, 0, 1], count=[2, 3, 1])

    sliced_ar = ar[::]
    assert sliced_ar.GetDimensionCount() == ar.GetDimensionCount()
    assert sliced_ar.Read() == ar.Read()

    sliced_ar = ar[0::]
    assert sliced_ar.GetDimensionCount() == ar.GetDimensionCount()
    assert sliced_ar.Read() == ar.Read()

    sliced_ar = ar[0:2:1]
    assert sliced_ar.GetDimensionCount() == ar.GetDimensionCount()
    assert sliced_ar.Read() == ar.Read()

    sliced_ar = ar[0:-1:1]
    assert sliced_ar.GetDimensionCount() == ar.GetDimensionCount()
    assert sliced_ar.Read() == ar.Read(array_start_idx=[0, 0, 0], count=[1, 3, 4])

    sliced_ar = ar[1:2:1]
    assert sliced_ar
    assert sliced_ar.GetDimensionCount() == ar.GetDimensionCount()
    assert sliced_ar.Read() == ar.Read(array_start_idx=[1, 0, 0], count=[1, 3, 4])

    sliced_ar = ar[::-1, 0, 0]
    assert sliced_ar
    assert sliced_ar.GetDimensionCount() == 1
    assert sliced_ar.GetDimensions()[0].GetSize() == 2
    read = sliced_ar.Read()
    assert struct.unpack("B" * len(read), read) == (0, 12)[::-1]

    sliced_ar = ar[0, 0, ::-1]
    assert sliced_ar
    read = sliced_ar.Read()
    assert struct.unpack("B" * len(read), read) == (0, 1, 2, 3)[::-1]

    sliced_ar = ar[0, 0, 2::-1]
    assert sliced_ar
    read = sliced_ar.Read()
    assert struct.unpack("B" * len(read), read) == (0, 1, 2, 3)[2::-1]

    sliced_ar = ar[0, 0, 2:0:-1]
    assert sliced_ar
    read = sliced_ar.Read()
    assert struct.unpack("B" * len(read), read) == (0, 1, 2, 3)[2:0:-1]

    sliced_ar = ar[0, 0, 0:4:1]
    assert sliced_ar
    read = sliced_ar.Read()
    assert struct.unpack("B" * len(read), read) == (0, 1, 2, 3)[0:4:1]

    sliced_ar = ar[0, 0, 0:4:2]
    assert sliced_ar
    read = sliced_ar.Read()
    assert struct.unpack("B" * len(read), read) == (0, 1, 2, 3)[0:4:2]

    sliced_ar = ar[0, 0, 0:3:2]
    assert sliced_ar
    read = sliced_ar.Read()
    assert struct.unpack("B" * len(read), read) == (0, 1, 2, 3)[0:3:2]

    sliced_ar = ar[0, 0, 3:0:-2]
    assert sliced_ar
    read = sliced_ar.Read()
    assert struct.unpack("B" * len(read), read) == (0, 1, 2, 3)[3:0:-2]

    sliced_ar = ar[0, 0, 3:1:-2]
    assert sliced_ar
    read = sliced_ar.Read()
    assert struct.unpack("B" * len(read), read) == (0, 1, 2, 3)[3:1:-2]


def test_mem_md_band_as_mdarray():

    drv = gdal.GetDriverByName("MEM")

    def get_array():
        ds = drv.Create("my_ds", 10, 5, 2, gdal.GDT_UInt16)
        ds.SetGeoTransform([2, 0.1, 0, 49, 0, -0.1])
        sr = osr.SpatialReference()
        sr.ImportFromEPSG(32631)
        ds.SetSpatialRef(sr)
        band = ds.GetRasterBand(1)
        band.SetUnitType("foo")
        band.SetNoDataValue(2)
        band.SetOffset(1.5)
        band.SetScale(2.5)
        band.SetMetadataItem("FOO", "BAR")
        band.WriteRaster(0, 0, 10, 5, struct.pack("H", 1), 1, 1)
        band.WriteRaster(0, 1, 1, 1, struct.pack("H", 2))
        band.WriteRaster(1, 0, 1, 1, struct.pack("H", 3))
        return (band.AsMDArray(), band.ReadRaster())

    ar, expected_data = get_array()
    assert ar
    assert ar.GetView("[...]")
    assert ar.GetBlockSize() == [1, 10]
    assert ar.GetDimensionCount() == 2
    dims = ar.GetDimensions()
    assert dims[0].GetName() == "Y"
    assert dims[0].GetSize() == 5
    var_y = dims[0].GetIndexingVariable()
    assert var_y
    assert var_y.GetDimensions()[0].GetSize() == 5
    assert struct.unpack("d" * 5, var_y.Read()) == (48.95, 48.85, 48.75, 48.65, 48.55)
    assert dims[1].GetName() == "X"
    assert dims[1].GetSize() == 10
    var_x = dims[1].GetIndexingVariable()
    assert var_x
    assert var_x.GetDimensions()[0].GetSize() == 10
    assert struct.unpack("d" * 10, var_x.Read()) == (
        2.05,
        2.15,
        2.25,
        2.35,
        2.45,
        2.55,
        2.65,
        2.75,
        2.85,
        2.95,
    )
    assert ar.GetDataType().GetClass() == gdal.GEDTC_NUMERIC
    assert ar.GetDataType().GetNumericDataType() == gdal.GDT_UInt16
    assert ar.Read() == expected_data
    assert struct.unpack(
        "H" * 4, ar.Read(array_start_idx=[1, 1], count=[2, 2], array_step=[-1, -1])
    ) == (1, 2, 3, 1)
    assert struct.unpack(
        "H" * 4,
        ar[1::-1, 1::-1].Read(array_start_idx=[0, 0], count=[2, 2], array_step=[1, 1]),
    ) == (1, 2, 3, 1)
    assert struct.unpack(
        "H" * 4, ar.Read(array_start_idx=[0, 0], count=[2, 2], array_step=[1, 1])
    ) == (1, 3, 2, 1)
    assert struct.unpack(
        "H" * 4,
        ar[1::-1, 1::-1].Read(
            array_start_idx=[1, 1], count=[2, 2], array_step=[-1, -1]
        ),
    ) == (1, 3, 2, 1)
    assert ar.Write(expected_data) == gdal.CE_None
    assert ar.Read() == expected_data
    assert ar.GetUnit() == "foo"
    assert ar.GetNoDataValueAsDouble() == 2
    assert ar.GetOffset() == 1.5
    assert ar.GetScale() == 2.5
    assert ar.GetSpatialRef() is not None
    assert len(ar.GetAttributes()) == 1
    attr = ar.GetAttribute("FOO")
    assert attr
    assert attr.Read() == "BAR"

    def get_array():
        ds = drv.Create("my_ds", 10, 5, 2)
        band = ds.GetRasterBand(1)
        return (band.AsMDArray(), band.ReadRaster())

    ar, expected_data = get_array()
    assert ar.GetSpatialRef() is None
    assert ar.GetUnit() == ""
    assert ar.GetNoDataValueAsRaw() is None


def test_mem_md_array_as_classic_dataset():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim_y = rg.CreateDimension("dim_y", None, None, 2)
    dim_x = rg.CreateDimension("dim_x", None, None, 3)

    dim_x_var = rg.CreateMDArray(
        "dim_x", [dim_x], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    dim_x_var.Write(struct.pack("d" * 3, 1.25, 2.25, 3.25))
    dim_x.SetIndexingVariable(dim_x_var)

    dim_y_var = rg.CreateMDArray(
        "dim_y", [dim_y], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    dim_y_var.Write(struct.pack("d" * 2, 10, 8))
    dim_y.SetIndexingVariable(dim_y_var)

    ar = rg.CreateMDArray("nodim", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte))
    with gdal.quiet_errors():
        assert not ar.AsClassicDataset(0, 0)

    ar = rg.CreateMDArray("1d", [dim_x], gdal.ExtendedDataType.Create(gdal.GDT_Byte))
    with gdal.quiet_errors():
        assert not ar.AsClassicDataset(1, 0)
    ds = ar.AsClassicDataset(0, 0)
    assert ds.RasterXSize == 3
    assert ds.RasterYSize == 1
    assert ds.RasterCount == 1
    assert not ds.GetSpatialRef()
    data = struct.pack("B" * 3, 0, 1, 2)
    assert ar.Write(data) == gdal.CE_None
    band = ds.GetRasterBand(1)
    assert len(band.ReadRaster(buf_type=gdal.GDT_UInt16)) == len(data) * 2
    assert band.ReadRaster(buf_type=gdal.GDT_UInt16) == struct.pack("H" * 3, 0, 1, 2)
    assert (
        band.WriteRaster(
            0, 0, 3, 1, struct.pack("H" * 3, 0, 1, 2), buf_type=gdal.GDT_UInt16
        )
        == gdal.CE_None
    )
    assert band.ReadRaster() == data

    ar = rg.CreateMDArray(
        "2d_string", [dim_y, dim_x], gdal.ExtendedDataType.CreateString()
    )
    with gdal.quiet_errors():
        assert not ar.AsClassicDataset(0, 1)

    # 2D
    ar = rg.CreateMDArray(
        "2d", [dim_y, dim_x], gdal.ExtendedDataType.Create(gdal.GDT_UInt16)
    )
    attr = ar.CreateAttribute(
        "attr_float64", [1], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    attr.Write(1.25)
    attr = ar.CreateAttribute("attr_strings", [2], gdal.ExtendedDataType.CreateString())
    attr.Write(["foo", "bar"])
    with gdal.quiet_errors():
        assert not ar.AsClassicDataset(0, 0)
        assert not ar.AsClassicDataset(0, 2)
        assert not ar.AsClassicDataset(2, 0)
    ds = ar.AsClassicDataset(1, 0)
    assert ds
    assert ds.RasterXSize == 3
    assert ds.RasterYSize == 2
    assert ds.RasterCount == 1
    assert ds.GetGeoTransform() == (0.75, 1.0, 0.0, 11.0, 0.0, -2.0)
    assert ds.GetMetadata() == {"attr_float64": "1.25", "attr_strings": "{foo,bar}"}
    band = ds.GetRasterBand(1)
    assert band.DataType == gdal.GDT_UInt16
    assert band.GetBlockSize() == [3, 1]
    data = struct.pack("H" * 6, 0, 1, 2, 3, 4, 5)
    assert ar.Write(data) == gdal.CE_None
    assert len(band.ReadRaster()) == len(data)
    assert band.ReadRaster() == data
    assert band.WriteRaster(0, 0, 3, 2, data) == gdal.CE_None
    assert band.ReadRaster() == data
    assert band.ReadBlock(0, 0) == struct.pack("H" * 3, 0, 1, 2)
    assert band.ReadBlock(0, 1) == struct.pack("H" * 3, 3, 4, 5)
    assert not band.GetNoDataValue()
    assert not band.GetScale()
    assert not band.GetOffset()
    assert band.GetUnitType() == ""

    assert ds.GetSpatialRef() is None
    sr = osr.SpatialReference()
    sr.ImportFromEPSG(32631)
    sr.SetDataAxisToSRSAxisMapping([2, 1])
    ar.SetSpatialRef(sr)
    got_sr = ds.GetSpatialRef()
    assert got_sr is not None
    assert got_sr.GetDataAxisToSRSAxisMapping() == [1, 2]
    got_sr.SetDataAxisToSRSAxisMapping([2, 1])
    assert got_sr.IsSame(sr)
    ar.SetNoDataValueDouble(2)
    ar.SetScale(1.5)
    ar.SetOffset(2.5)
    ar.SetUnit("foo")
    assert band.GetNoDataValue() == 2
    assert band.GetScale() == 1.5
    assert band.GetOffset() == 2.5
    assert band.GetUnitType() == "foo"

    # 3D with band dimension first
    dim_bands = rg.CreateDimension("dim_bands", None, None, 2)
    dim_bands_var = rg.CreateMDArray(
        "dim_bands", [dim_bands], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    dim_bands_var.Write(struct.pack("d" * 2, 1.25, 2.25))
    dim_bands_var.SetUnit("my_unit")
    dim_bands.SetIndexingVariable(dim_bands_var)

    ar = rg.CreateMDArray(
        "3d_band_first",
        [dim_bands, dim_y, dim_x],
        gdal.ExtendedDataType.Create(gdal.GDT_UInt16),
    )
    ds = ar.AsClassicDataset(2, 1)
    assert ds
    assert ds.RasterXSize == 3
    assert ds.RasterYSize == 2
    assert ds.RasterCount == 2
    band = ds.GetRasterBand(2)
    assert band.GetMetadata() == {
        "DIM_dim_bands_INDEX": "1",
        "DIM_dim_bands_UNIT": "my_unit",
        "DIM_dim_bands_VALUE": "2.25",
    }
    assert band.GetBlockSize() == [3, 1]
    data = struct.pack("H" * 6, 0, 1, 2, 3, 4, 5)
    assert ar.Write(data, array_start_idx=[1, 0, 0], count=[1, 2, 3]) == gdal.CE_None
    assert band.ReadRaster() == data
    assert band.WriteRaster(0, 0, 3, 2, data) == gdal.CE_None
    assert band.ReadRaster() == data

    # 3D with band dimension last
    ar = rg.CreateMDArray(
        "3d_band_last",
        [dim_y, dim_x, dim_bands],
        gdal.ExtendedDataType.Create(gdal.GDT_UInt16),
    )
    ds = ar.AsClassicDataset(1, 0)
    assert ds
    assert ds.RasterXSize == 3
    assert ds.RasterYSize == 2
    assert ds.RasterCount == 2
    band = ds.GetRasterBand(2)
    assert band.GetMetadata() == {
        "DIM_dim_bands_INDEX": "1",
        "DIM_dim_bands_UNIT": "my_unit",
        "DIM_dim_bands_VALUE": "2.25",
    }
    assert band.GetBlockSize() == [3, 1]
    data = struct.pack("H" * 6, 0, 1, 2, 3, 4, 5)
    assert ar.Write(data, array_start_idx=[0, 0, 1], count=[2, 3, 1]) == gdal.CE_None
    assert band.ReadRaster() == data
    assert band.WriteRaster(0, 0, 3, 2, data) == gdal.CE_None
    assert band.ReadRaster() == data

    # 4D
    dim_time_subset = rg.CreateDimension("subset_time_1_2_3", None, None, 3)
    dim_time_subset_var = rg.CreateMDArray(
        "subset_time_1_2_3",
        [dim_time_subset],
        gdal.ExtendedDataType.Create(gdal.GDT_Float64),
    )
    dim_time_subset_var.Write(struct.pack("d" * 3, 1.5, 2.5, 3.5))
    dim_time_subset.SetIndexingVariable(dim_time_subset_var)
    ar = rg.CreateMDArray(
        "4d",
        [dim_time_subset, dim_y, dim_x, dim_bands],
        gdal.ExtendedDataType.Create(gdal.GDT_UInt16),
    )
    ds = ar.AsClassicDataset(2, 1)
    band = ds.GetRasterBand(5)
    assert band.GetMetadata() == {
        "DIM_dim_bands_INDEX": "0",
        "DIM_dim_bands_UNIT": "my_unit",
        "DIM_dim_bands_VALUE": "1.25",
        "DIM_time_INDEX": "5",
        "DIM_time_VALUE": "3.5",
    }


def test_mem_md_array_transpose():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim_z = rg.CreateDimension("dim_z", None, None, 4)
    dim_y = rg.CreateDimension("dim_y", None, None, 2)
    dim_x = rg.CreateDimension("dim_x", None, None, 3)
    ar = rg.CreateMDArray(
        "ar", [dim_z, dim_y, dim_x], gdal.ExtendedDataType.Create(gdal.GDT_UInt16)
    )

    attr = ar.CreateAttribute(
        "attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    assert attr.Write(1) == gdal.CE_None

    data = array.array("H", list(range(24))).tobytes()
    assert ar.Write(data) == gdal.CE_None

    with gdal.quiet_errors():
        assert not ar.Transpose([])  # 0 axis
        assert not ar.Transpose([0, 1])  # missing axis
        assert not ar.Transpose([0, 1, 2, 3])  # too many axis
        assert not ar.Transpose([0, 1, 3])  # invalid axis number
        assert not ar.Transpose([0, 1, -2])  # invalid axis number
        assert not ar.Transpose([0, 1, -1])  # missing axis
        assert not ar.Transpose([0, 1, 1])  # repeated axis

    # Idendity
    transposed = ar.Transpose([0, 1, 2])
    assert transposed.GetDimensionCount() == 3
    dims = transposed.GetDimensions()
    assert dims[0].GetName() == "dim_z"
    assert dims[1].GetName() == "dim_y"
    assert dims[2].GetName() == "dim_x"
    assert transposed.Read() == data

    assert transposed.GetUnit() == ar.GetUnit()
    assert transposed.GetScale() == ar.GetScale()
    assert transposed.GetOffset() == ar.GetOffset()
    assert transposed.GetNoDataValueAsRaw() == ar.GetNoDataValueAsRaw()
    assert transposed.GetSpatialRef() is None
    assert transposed.GetAttribute("attr") is not None
    assert len(transposed.GetAttributes()) == 1

    assert ar.SetUnit("foo") == gdal.CE_None
    assert ar.SetScale(1) == gdal.CE_None
    assert ar.SetOffset(2) == gdal.CE_None
    assert ar.SetNoDataValueDouble(1) == gdal.CE_None
    sr = osr.SpatialReference()
    sr.ImportFromEPSG(32631)
    sr.SetDataAxisToSRSAxisMapping([3, 2])
    assert ar.SetSpatialRef(sr) == gdal.CE_None
    assert ar.GetSpatialRef()

    assert transposed.GetUnit() == ar.GetUnit()
    assert transposed.GetScale() == ar.GetScale()
    assert transposed.GetOffset() == ar.GetOffset()
    assert transposed.GetNoDataValueAsRaw() == ar.GetNoDataValueAsRaw()
    assert transposed.GetSpatialRef().IsSame(ar.GetSpatialRef())

    # Idendity with one extra axis
    transposed = ar.Transpose([0, 1, -1, 2])
    assert transposed.GetDimensionCount() == 4
    dims = transposed.GetDimensions()
    assert dims[0].GetName() == "dim_z"
    assert dims[1].GetName() == "dim_y"
    assert dims[2].GetName() == "newaxis"
    assert dims[3].GetName() == "dim_x"
    assert transposed.Read() == data
    assert transposed.GetBlockSize() == [0, 0, 0, 0]

    got_sr = transposed.GetSpatialRef()
    assert got_sr
    assert got_sr.GetDataAxisToSRSAxisMapping() == [4, 2]

    # Full transpose
    transposed = ar.Transpose([2, 1, 0])
    assert transposed.GetDimensionCount() == 3
    dims = transposed.GetDimensions()
    assert dims[0].GetName() == "dim_x"
    assert dims[1].GetName() == "dim_y"
    assert dims[2].GetName() == "dim_z"
    transposed_data = transposed.Read()
    assert struct.unpack("H" * 24, transposed_data) == (
        0,
        6,
        12,
        18,
        3,
        9,
        15,
        21,
        1,
        7,
        13,
        19,
        4,
        10,
        16,
        22,
        2,
        8,
        14,
        20,
        5,
        11,
        17,
        23,
    )
    assert transposed.Read(buffer_stride=[1, 3, 6]) == data
    assert transposed.Write(struct.pack("H", 0) * 24) == gdal.CE_None
    assert ar.Read() == struct.pack("H", 0) * 24
    assert transposed.Write(transposed_data) == gdal.CE_None
    assert ar.Read() == data

    # Rotation of axis
    transposed = ar.Transpose([1, 2, 0])
    assert transposed.Read(buffer_stride=[3, 1, 6]) == data


def test_mem_md_array_single_dim_non_contiguous_copy():

    drv = gdal.GetDriverByName("MEM")
    nvalues = 30
    spacing = 63
    data = array.array("B", list(range(nvalues))).tobytes()
    for t in (
        gdal.GDT_Byte,
        gdal.GDT_Int16,
        gdal.GDT_Int32,
        gdal.GDT_Float64,
        gdal.GDT_CFloat64,
    ):
        ds = drv.CreateMultiDimensional("myds")
        rg = ds.GetRootGroup()
        dim = rg.CreateDimension("dim0", "", "", nvalues)
        ar = rg.CreateMDArray("ar", [dim], gdal.ExtendedDataType.Create(t))
        assert ar
        assert (
            ar.Write(data, buffer_datatype=gdal.ExtendedDataType.Create(gdal.GDT_Byte))
            == gdal.CE_None
        )
        got_data = ar.Read(
            buffer_stride=[spacing],
            buffer_datatype=gdal.ExtendedDataType.Create(gdal.GDT_Byte),
        )
        assert len(got_data) == (nvalues - 1) * spacing + 1
        got_data = struct.unpack("B" * len(got_data), got_data)
        for i in range(nvalues):
            assert got_data[i * spacing] == i


def test_mem_md_array_get_unscaled_0dim():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    myarray = rg.CreateMDArray(
        "myarray", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    assert myarray

    data = array.array("B", [1]).tobytes()
    assert myarray.Write(data) == gdal.CE_None

    myarray.SetOffset(1.5)
    myarray.SetScale(200.5)

    unscaled = myarray.GetUnscaled()
    assert unscaled.GetDataType().GetNumericDataType() == gdal.GDT_Float64
    assert struct.unpack("d" * 1, unscaled.Read())[0] == 1 * 200.5 + 1.5

    float32dt = gdal.ExtendedDataType.Create(gdal.GDT_Float32)
    assert (
        struct.unpack("f" * 1, unscaled.Read(buffer_datatype=float32dt))[0]
        == 1 * 200.5 + 1.5
    )

    assert unscaled.Write(struct.pack("d" * 1, 2 * 200.5 + 1.5)) == gdal.CE_None
    assert struct.unpack("B" * 1, myarray.Read())[0] == 2

    assert unscaled.Write(struct.pack("d" * 1, 2.1 * 200.5 + 1.5)) == gdal.CE_None
    assert struct.unpack("B" * 1, myarray.Read())[0] == 2

    assert unscaled.Write(struct.pack("d" * 1, 1.9 * 200.5 + 1.5)) == gdal.CE_None
    assert struct.unpack("B" * 1, myarray.Read())[0] == 2

    assert (
        unscaled.Write(struct.pack("f" * 1, 3 * 200.5 + 1.5), buffer_datatype=float32dt)
        == gdal.CE_None
    )
    assert struct.unpack("B" * 1, myarray.Read())[0] == 3


def test_mem_md_array_get_unscaled_0dim_complex():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    myarray = rg.CreateMDArray(
        "myarray", [], gdal.ExtendedDataType.Create(gdal.GDT_CInt16)
    )
    assert myarray

    data = array.array("H", [1, 2]).tobytes()
    assert myarray.Write(data) == gdal.CE_None

    myarray.SetOffset(1.5)
    myarray.SetScale(200.5)

    unscaled = myarray.GetUnscaled()
    assert unscaled.GetDataType().GetNumericDataType() == gdal.GDT_CFloat64
    assert struct.unpack("d" * 2, unscaled.Read()) == (1 * 200.5 + 1.5, 2 * 200.5 + 1.5)

    cfloat32dt = gdal.ExtendedDataType.Create(gdal.GDT_CFloat32)
    assert struct.unpack("f" * 2, unscaled.Read(buffer_datatype=cfloat32dt)) == (
        1 * 200.5 + 1.5,
        2 * 200.5 + 1.5,
    )

    assert (
        unscaled.Write(struct.pack("d" * 2, 3 * 200.5 + 1.5, 4 * 200.5 + 1.5))
        == gdal.CE_None
    )
    assert struct.unpack("H" * 2, myarray.Read()) == (3, 4)

    assert (
        unscaled.Write(
            struct.pack("f" * 2, 5 * 200.5 + 1.5, 6 * 200.5 + 1.5),
            buffer_datatype=cfloat32dt,
        )
        == gdal.CE_None
    )
    assert struct.unpack("H" * 2, myarray.Read()) == (5, 6)


def test_mem_md_array_get_unscaled_0dim_non_matching_nodata():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    myarray = rg.CreateMDArray(
        "myarray", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    assert myarray

    data = array.array("B", [1]).tobytes()
    assert myarray.Write(data) == gdal.CE_None

    myarray.SetOffset(1.5)
    myarray.SetScale(200.5)
    myarray.SetNoDataValueDouble(3)

    unscaled = myarray.GetUnscaled()
    assert unscaled.GetDataType().GetNumericDataType() == gdal.GDT_Float64
    nodata = unscaled.GetNoDataValueAsDouble()
    assert math.isnan(nodata)
    assert struct.unpack("d" * 1, unscaled.Read())[0] == 1 * 200.5 + 1.5
    assert (
        struct.unpack(
            "f" * 1,
            unscaled.Read(
                buffer_datatype=gdal.ExtendedDataType.Create(gdal.GDT_Float32)
            ),
        )[0]
        == 1 * 200.5 + 1.5
    )

    assert unscaled.Write(struct.pack("d" * 1, 2 * 200.5 + 1.5)) == gdal.CE_None
    assert struct.unpack("B" * 1, myarray.Read())[0] == 2


def test_mem_md_array_get_unscaled_0dim_matching_nodata():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    myarray = rg.CreateMDArray(
        "myarray", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    assert myarray

    data = array.array("B", [1]).tobytes()
    assert myarray.Write(data) == gdal.CE_None

    myarray.SetOffset(1.5)
    myarray.SetScale(200.5)
    myarray.SetNoDataValueDouble(1)

    unscaled = myarray.GetUnscaled()
    assert unscaled.GetDataType().GetNumericDataType() == gdal.GDT_Float64
    nodata = unscaled.GetNoDataValueAsDouble()
    assert math.isnan(nodata)
    assert math.isnan(struct.unpack("d" * 1, unscaled.Read())[0])

    assert unscaled.Write(struct.pack("d" * 1, 2 * 200.5 + 1.5)) == gdal.CE_None
    assert struct.unpack("B" * 1, myarray.Read())[0] == 2

    assert unscaled.Write(struct.pack("d" * 1, nodata)) == gdal.CE_None
    assert struct.unpack("B" * 1, myarray.Read())[0] == 1


def test_mem_md_array_get_unscaled_0dim_matching_nodata_complex():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    myarray = rg.CreateMDArray(
        "myarray", [], gdal.ExtendedDataType.Create(gdal.GDT_CInt16)
    )
    assert myarray

    data = array.array("H", [1, 2]).tobytes()
    assert myarray.Write(data) == gdal.CE_None

    myarray.SetOffset(1.5)
    myarray.SetScale(200.5)
    myarray.SetNoDataValueRaw(struct.pack("H" * 2, 1, 2))

    unscaled = myarray.GetUnscaled()
    assert unscaled.GetDataType().GetNumericDataType() == gdal.GDT_CFloat64
    nodata = unscaled.GetNoDataValueAsDouble()
    assert math.isnan(nodata)
    assert math.isnan(struct.unpack("d" * 2, unscaled.Read())[0])
    assert math.isnan(struct.unpack("d" * 2, unscaled.Read())[1])

    assert unscaled.Write(struct.pack("d" * 2, nodata, nodata)) == gdal.CE_None
    assert struct.unpack("H" * 2, myarray.Read()) == (1, 2)


def test_mem_md_array_get_unscaled_3dim():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim0 = rg.CreateDimension("dim0", None, None, 2)
    dim1 = rg.CreateDimension("dim1", None, None, 3)
    dim2 = rg.CreateDimension("dim2", None, None, 4)
    myarray = rg.CreateMDArray(
        "myarray", [dim0, dim1, dim2], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    assert myarray

    data = array.array("B", list(range(24))).tobytes()
    assert myarray.Write(data) == gdal.CE_None

    assert myarray.GetUnscaled().Read() == myarray.Read()

    myarray.SetOffset(1.5)
    myarray.SetScale(200.5)

    unscaled = myarray.GetUnscaled()
    assert unscaled.GetOffset() is None
    assert unscaled.GetScale() is None
    assert unscaled.GetNoDataValueAsRaw() is None
    assert unscaled.GetSpatialRef() is None
    assert unscaled.GetUnit() == myarray.GetUnit()
    assert unscaled.GetBlockSize() == myarray.GetBlockSize()
    assert [x.GetSize() for x in unscaled.GetDimensions()] == [
        x.GetSize() for x in myarray.GetDimensions()
    ]
    assert unscaled.GetDataType().GetNumericDataType() == gdal.GDT_Float64
    expected_data = [x * 200.5 + 1.5 for x in range(24)]
    unscaled_data = unscaled.Read()
    assert [x for x in struct.unpack("d" * 24, unscaled_data)] == expected_data
    float32_dt = gdal.ExtendedDataType.Create(gdal.GDT_Float32)
    unscaled_data_float32 = unscaled.Read(buffer_datatype=float32_dt)
    assert [x for x in struct.unpack("f" * 24, unscaled_data_float32)] == expected_data

    assert myarray.Write(b"\x00" * 24) == gdal.CE_None
    assert myarray.Read() != data

    assert unscaled.Write(unscaled_data) == gdal.CE_None
    assert myarray.Read() == data

    assert myarray.Write(b"\x00" * 24) == gdal.CE_None
    assert myarray.Read() != data

    assert (
        unscaled.Write(unscaled_data_float32, buffer_datatype=float32_dt)
        == gdal.CE_None
    )
    assert myarray.Read() == data

    myarray.SetNoDataValueDouble(1)
    unscaled = myarray.GetUnscaled()
    assert math.isnan(unscaled.GetNoDataValueAsDouble())
    unscaled_data_with_nan = unscaled.Read()
    got_data = [x for x in struct.unpack("d" * 24, unscaled_data_with_nan)]
    expected_data = [
        float("nan") if x == 1 else x * 200.5 + 1.5
        for x in struct.unpack("B" * 24, myarray.Read())
    ]
    for i in range(24):
        if math.isnan(expected_data[i]):
            assert math.isnan(got_data[i])
        else:
            assert got_data[i] == expected_data[i]

    assert myarray.Write(b"\x00" * 24) == gdal.CE_None
    assert myarray.Read() != data

    assert unscaled.Write(unscaled_data_with_nan) == gdal.CE_None
    assert myarray.Read() == data


def test_mem_md_array_get_unscaled_1dim_complex():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim0 = rg.CreateDimension("dim0", None, None, 2)
    myarray = rg.CreateMDArray(
        "myarray", [dim0], gdal.ExtendedDataType.Create(gdal.GDT_CInt16)
    )
    assert myarray

    data = array.array("H", [1, 2, 3, 4]).tobytes()
    assert myarray.Write(data) == gdal.CE_None

    assert myarray.GetUnscaled().Read() == myarray.Read()

    assert myarray.GetUnscaled().AdviseRead() == gdal.CE_None

    myarray.SetOffset(1.5)
    myarray.SetScale(200.5)
    myarray.SetNoDataValueRaw(struct.pack("H" * 2, 1, 2))

    unscaled = myarray.GetUnscaled()
    assert unscaled.GetDataType().GetNumericDataType() == gdal.GDT_CFloat64
    unscaled_data_with_nan = unscaled.Read()
    got_data = [x for x in struct.unpack("d" * 4, unscaled_data_with_nan)]
    assert math.isnan(got_data[0])
    assert math.isnan(got_data[1])
    assert got_data[2] == 3 * 200.5 + 1.5
    assert got_data[3] == 4 * 200.5 + 1.5

    assert unscaled.Write(unscaled_data_with_nan) == gdal.CE_None
    assert myarray.Read() == data


def test_mem_md_array_get_mask():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    myarray = rg.CreateMDArray(
        "myarray_emptydim", [], gdal.ExtendedDataType.Create(gdal.GDT_Int16)
    )
    mask = myarray.GetMask()
    assert mask is not None
    assert struct.unpack("B", mask.Read())[0] == 1
    assert (
        struct.unpack(
            "H", mask.Read(buffer_datatype=gdal.ExtendedDataType.Create(gdal.GDT_Int16))
        )[0]
        == 1
    )

    myarray.SetNoDataValueDouble(0)
    assert struct.unpack("B", mask.Read())[0] == 0
    assert (
        struct.unpack(
            "H", mask.Read(buffer_datatype=gdal.ExtendedDataType.Create(gdal.GDT_Int16))
        )[0]
        == 0
    )

    assert myarray.GetMask().AdviseRead() == gdal.CE_None

    dim0 = rg.CreateDimension("dim0", None, None, 2)
    dim1 = rg.CreateDimension("dim1", None, None, 3)
    dim2 = rg.CreateDimension("dim2", None, None, 4)
    myarray = rg.CreateMDArray(
        "myarray_string", [dim0], gdal.ExtendedDataType.CreateString()
    )
    # Non-numeric array unsupported
    with gdal.quiet_errors():
        assert not myarray.GetMask()

    myarray = rg.CreateMDArray(
        "myarray", [dim0, dim1, dim2], gdal.ExtendedDataType.Create(gdal.GDT_Int32)
    )
    data = array.array("I", list(range(24))).tobytes()
    assert myarray.Write(data) == gdal.CE_None

    mask = myarray.GetMask()
    assert mask.GetOffset() is None
    assert mask.GetScale() is None
    with gdal.quiet_errors():
        assert mask.Write(mask.Read()) == gdal.CE_Failure
    assert mask.GetNoDataValueAsRaw() is None
    assert mask.GetSpatialRef() is None
    assert mask.GetUnit() == myarray.GetUnit()
    assert mask.GetBlockSize() == myarray.GetBlockSize()
    assert [x.GetSize() for x in mask.GetDimensions()] == [
        x.GetSize() for x in myarray.GetDimensions()
    ]
    assert mask.GetDataType().GetNumericDataType() == gdal.GDT_Byte
    # Case when we don't need to read the underlying array at all: the mask is always valid
    assert [x for x in struct.unpack("B" * 24, mask.Read())] == [1] * 24
    assert [x for x in struct.unpack("B" * 24, mask.Read(buffer_stride=[1, 2, 6]))] == [
        1
    ] * 24
    assert [
        x
        for x in struct.unpack(
            "H" * 24,
            mask.Read(buffer_datatype=gdal.ExtendedDataType.Create(gdal.GDT_Int16)),
        )
    ] == [1] * 24
    assert [
        x
        for x in struct.unpack(
            "H" * 24,
            mask.Read(
                buffer_datatype=gdal.ExtendedDataType.Create(gdal.GDT_Int16),
                buffer_stride=[1, 2, 6],
            ),
        )
    ] == [1] * 24

    # Test no data value
    myarray.SetNoDataValueDouble(10)
    mask = myarray.GetMask()
    expected_data = [1] * 24
    expected_data[10] = 0
    assert [x for x in struct.unpack("B" * 24, mask.Read())] == expected_data
    assert [x for x in struct.unpack("B" * 24, mask.Read(buffer_stride=[1, 2, 6]))] == [
        1,
        1,
        1,
        1,
        1,
        1,
        1,
        1,
        1,
        1,
        1,
        1,
        1,
        1,
        1,
        1,
        0,
        1,
        1,
        1,
        1,
        1,
        1,
        1,
    ]
    assert [
        x
        for x in struct.unpack(
            "H" * 24,
            mask.Read(buffer_datatype=gdal.ExtendedDataType.Create(gdal.GDT_Int16)),
        )
    ] == expected_data

    # Test missing_value, _FillValue, valid_min, valid_max
    bytedt = gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    attr = myarray.CreateAttribute("missing_value", [1], bytedt)
    assert attr.Write(8) == gdal.CE_None
    attr = myarray.CreateAttribute("_FillValue", [1], bytedt)
    assert attr.Write(9) == gdal.CE_None
    attr = myarray.CreateAttribute("valid_min", [1], bytedt)
    assert attr.Write(2) == gdal.CE_None
    attr = myarray.CreateAttribute("valid_max", [1], bytedt)
    assert attr.Write(22) == gdal.CE_None
    mask = myarray.GetMask()
    expected_data = [1] * 24
    expected_data[0] = 0
    expected_data[1] = 0
    expected_data[8] = 0
    expected_data[9] = 0
    expected_data[10] = 0
    expected_data[23] = 0
    assert [x for x in struct.unpack("B" * 24, mask.Read())] == expected_data

    # Test valid_range
    myarray = rg.CreateMDArray(
        "myarray_valid_range",
        [dim0, dim1, dim2],
        gdal.ExtendedDataType.Create(gdal.GDT_Int16),
    )
    data = array.array("H", list(range(24))).tobytes()
    assert myarray.Write(data) == gdal.CE_None
    attr = myarray.CreateAttribute("valid_range", [2], bytedt)
    assert attr.Write([1, 22]) == gdal.CE_None
    mask = myarray.GetMask()
    expected_data = [1] * 24
    expected_data[0] = 0
    expected_data[23] = 0
    assert [x for x in struct.unpack("B" * 24, mask.Read())] == expected_data

    try:
        import numpy

        has_numpy = True
    except ImportError:
        has_numpy = False

    if has_numpy:
        ma = myarray.ReadAsMaskedArray()
        assert ma[0, 0, 0] is numpy.ma.masked
        assert ma[0, 0, 1] is not numpy.ma.masked

    # Test array with nan
    myarray = rg.CreateMDArray(
        "myarray_with_nan", [dim0], gdal.ExtendedDataType.Create(gdal.GDT_Float32)
    )
    assert myarray.Write(struct.pack("f" * 2, 0, float("nan"))) == gdal.CE_None

    mask = myarray.GetMask()
    assert [x for x in struct.unpack("B" * 2, mask.Read())] == [1, 0]

    # Test all data types
    for dt, v, nv, expected in [
        (gdal.GDT_Byte, 1, 1, [1, 0]),
        (gdal.GDT_Byte, 1, 1.5, [1, 1]),
        (gdal.GDT_Int16, 1, 1, [1, 0]),
        (gdal.GDT_UInt16, 1, 1, [1, 0]),
        (gdal.GDT_Int32, 1, 1, [1, 0]),
        (gdal.GDT_UInt32, 1, 1, [1, 0]),
        (gdal.GDT_Float32, 1, 1, [1, 0]),
        (gdal.GDT_Float32, 1.5, 1.5, [1, 0]),
        (gdal.GDT_Float64, 1, 1, [1, 0]),
        (gdal.GDT_Float64, 1.5, 1.5, [1, 0]),
        (gdal.GDT_CInt16, 1, 1, [1, 0]),
    ]:
        myarray = rg.CreateMDArray(
            "array_dt_" + gdal.GetDataTypeName(dt) + "_" + str(v) + "_" + str(nv),
            [dim0],
            gdal.ExtendedDataType.Create(dt),
        )
        assert (
            myarray.Write(
                struct.pack("d" * 2, 0, v),
                buffer_datatype=gdal.ExtendedDataType.Create(gdal.GDT_Float64),
            )
            == gdal.CE_None
        )
        myarray.SetNoDataValueDouble(nv)
        mask = myarray.GetMask()
        assert [
            x for x in struct.unpack("B" * 2, mask.Read())
        ] == expected, myarray.GetName()


@gdaltest.enable_exceptions()
def test_mem_md_array_get_mask_unmask_flags_option_flag_values_only():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    dim0 = rg.CreateDimension("dim0", None, None, 2)
    dim1 = rg.CreateDimension("dim1", None, None, 3)
    myarray = rg.CreateMDArray(
        "myarray", [dim0, dim1], gdal.ExtendedDataType.Create(gdal.GDT_Int16)
    )

    with pytest.raises(Exception, match="no flag_meanings attribute"):
        myarray.GetMask(["UNMASK_FLAGS=bar"])

    flag_meanings = myarray.CreateAttribute(
        "flag_meanings", [], gdal.ExtendedDataType.CreateString()
    )

    with pytest.raises(Exception, match="Cannot read flag_meanings attribute"):
        myarray.GetMask(["UNMASK_FLAGS=bar"])

    flag_meanings.Write("one two _255")
    with pytest.raises(
        Exception, match="Cannot find flag_values and/or flag_masks attribute"
    ):
        myarray.GetMask(["UNMASK_FLAGS=one"])

    flag_values = myarray.CreateAttribute(
        "flag_values", [1], gdal.ExtendedDataType.Create(gdal.GDT_Int16)
    )
    with pytest.raises(
        Exception,
        match="Number of values in flag_values attribute is different from the one in flag_meanings",
    ):
        myarray.GetMask(["UNMASK_FLAGS=one"])

    myarray.DeleteAttribute("flag_values")
    flag_values = myarray.CreateAttribute(
        "flag_values", [3], gdal.ExtendedDataType.Create(gdal.GDT_Int16)
    )
    flag_values.Write(array.array("h", [1, 2, 255]))

    with pytest.raises(Exception, match="Cannot fing flag three in flag_meanings"):
        myarray.GetMask(["UNMASK_FLAGS=three"])

    data = array.array("h", list(range(2 * 3))).tobytes()
    assert myarray.Write(data) == gdal.CE_None

    mask = myarray.GetMask(["UNMASK_FLAGS=two"])
    assert array.array("B", mask.Read()).tolist() == [0, 0, 1, 0, 0, 0]
    assert array.array(
        "h", mask.Read(buffer_datatype=gdal.ExtendedDataType.Create(gdal.GDT_Int16))
    ).tolist() == [0, 0, 1, 0, 0, 0]

    mask = myarray.GetMask(["UNMASK_FLAGS=two,one"])
    assert array.array("B", mask.Read()).tolist() == [0, 1, 1, 0, 0, 0]

    mask = myarray.GetMask(["UNMASK_FLAGS=_255"])
    assert array.array("B", mask.Read()).tolist() == [0, 0, 0, 0, 0, 0]


@gdaltest.enable_exceptions()
def test_mem_md_array_get_mask_unmask_flags_option_flag_masks_only():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    dim0 = rg.CreateDimension("dim0", None, None, 2)
    dim1 = rg.CreateDimension("dim1", None, None, 3)
    myarray = rg.CreateMDArray(
        "myarray", [dim0, dim1], gdal.ExtendedDataType.Create(gdal.GDT_Int16)
    )

    flag_meanings = myarray.CreateAttribute(
        "flag_meanings", [], gdal.ExtendedDataType.CreateString()
    )
    flag_meanings.Write("bit_0 bit_1 bit_7")

    flag_masks = myarray.CreateAttribute(
        "flag_masks", [1], gdal.ExtendedDataType.Create(gdal.GDT_Int16)
    )
    with pytest.raises(
        Exception,
        match="Number of values in flag_masks attribute is different from the one in flag_meanings",
    ):
        myarray.GetMask(["UNMASK_FLAGS=one"])

    myarray.DeleteAttribute("flag_masks")
    flag_masks = myarray.CreateAttribute(
        "flag_masks", [3], gdal.ExtendedDataType.Create(gdal.GDT_Int16)
    )
    flag_masks.Write(array.array("h", [1 << 0, 1 << 1, 1 << 7]))

    data = array.array("h", list(range(2 * 3))).tobytes()
    assert myarray.Write(data) == gdal.CE_None

    mask = myarray.GetMask(["UNMASK_FLAGS=bit_0"])
    assert array.array("B", mask.Read()).tolist() == [0, 1, 0, 1, 0, 1]

    mask = myarray.GetMask(["UNMASK_FLAGS=bit_1"])
    assert array.array("B", mask.Read()).tolist() == [0, 0, 1, 1, 0, 0]

    mask = myarray.GetMask(["UNMASK_FLAGS=bit_1,bit_0"])
    assert array.array("B", mask.Read()).tolist() == [0, 1, 1, 1, 0, 1]


@gdaltest.enable_exceptions()
def test_mem_md_array_get_mask_unmask_flags_option_flag_values_and_masks():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    dim0 = rg.CreateDimension("dim0", None, None, 2)
    dim1 = rg.CreateDimension("dim1", None, None, 3)
    myarray = rg.CreateMDArray(
        "myarray", [dim0, dim1], gdal.ExtendedDataType.Create(gdal.GDT_Int16)
    )

    flag_meanings = myarray.CreateAttribute(
        "flag_meanings", [], gdal.ExtendedDataType.CreateString()
    )
    flag_meanings.Write("valid invalid_1 invalid_2 invalid_3")

    flag_values = myarray.CreateAttribute(
        "flag_values", [4], gdal.ExtendedDataType.Create(gdal.GDT_Int16)
    )
    flag_values.Write(array.array("h", [1, (1 << 1), (2 << 1), (3 << 1)]))

    flag_masks = myarray.CreateAttribute(
        "flag_masks", [4], gdal.ExtendedDataType.Create(gdal.GDT_Int16)
    )
    flag_masks.Write(array.array("h", [1, 6, 6, 6]))

    data = array.array("h", [1, (1 << 1), (2 << 1), (3 << 1), 1, 1]).tobytes()
    assert myarray.Write(data) == gdal.CE_None

    mask = myarray.GetMask(["UNMASK_FLAGS=valid"])
    assert array.array("B", mask.Read()).tolist() == [1, 0, 0, 0, 1, 1]

    mask = myarray.GetMask(["UNMASK_FLAGS=invalid_1"])
    assert array.array("B", mask.Read()).tolist() == [0, 1, 0, 0, 0, 0]


def test_mem_md_array_resolvemdarray():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    a = rg.CreateGroup("a")
    aa = a.CreateGroup("aa")
    a.CreateGroup("ab")
    b = rg.CreateGroup("b")

    a.CreateMDArray("var_a", [], gdal.ExtendedDataType.Create(gdal.GDT_Int16))
    a.CreateMDArray("var_c", [], gdal.ExtendedDataType.Create(gdal.GDT_Int16))
    aa.CreateMDArray("var_aa", [], gdal.ExtendedDataType.Create(gdal.GDT_Int16))
    b.CreateMDArray("var_b", [], gdal.ExtendedDataType.Create(gdal.GDT_Int16))
    b.CreateMDArray("var_c", [], gdal.ExtendedDataType.Create(gdal.GDT_Int16))

    assert rg.ResolveMDArray("x", "/") is None
    with gdal.ExceptionMgr(), pytest.raises(Exception):
        rg.ResolveMDArray("x", "/")

    assert rg.ResolveMDArray("/a/var_a", "/").GetFullName() == "/a/var_a"
    assert rg.ResolveMDArray("var_a", "/").GetFullName() == "/a/var_a"
    assert rg.ResolveMDArray("var_a", "").GetFullName() == "/a/var_a"
    assert rg.ResolveMDArray("var_a", "/a").GetFullName() == "/a/var_a"
    assert rg.ResolveMDArray("var_a", "/a/aa").GetFullName() == "/a/var_a"
    assert rg.ResolveMDArray("var_a", "/a/ab").GetFullName() == "/a/var_a"
    assert rg.ResolveMDArray("var_a", "/b").GetFullName() == "/a/var_a"

    assert rg.ResolveMDArray("var_aa", "/").GetFullName() == "/a/aa/var_aa"
    assert a.ResolveMDArray("var_aa", "/").GetFullName() == "/a/aa/var_aa"
    assert a.ResolveMDArray("var_aa", "/aa").GetFullName() == "/a/aa/var_aa"
    assert aa.ResolveMDArray("var_aa", "/").GetFullName() == "/a/aa/var_aa"

    assert rg.ResolveMDArray("var_b", "").GetFullName() == "/b/var_b"

    assert rg.ResolveMDArray("var_c", "/a").GetFullName() == "/a/var_c"
    assert rg.ResolveMDArray("var_c", "/a/aa").GetFullName() == "/a/var_c"
    assert rg.ResolveMDArray("var_c", "/b").GetFullName() == "/b/var_c"
    assert a.ResolveMDArray("var_c", "").GetFullName() == "/a/var_c"
    assert b.ResolveMDArray("var_c", "").GetFullName() == "/b/var_c"


def test_mem_md_array_statistics():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim0 = rg.CreateDimension("dim0", "unspecified type", "unspecified direction", 2)
    dim1 = rg.CreateDimension("dim1", "unspecified type", "unspecified direction", 3)
    float64dt = gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    ar = rg.CreateMDArray("myarray", [dim0, dim1], float64dt)
    ar.SetNoDataValueDouble(6)
    data = struct.pack("d" * 6, 1, 2, 3, 4, 5, 6)
    ar.Write(data)

    stats = ar.ComputeStatistics(False)
    assert stats.min == 1.0
    assert stats.max == 5.0
    assert stats.mean == 3.0
    assert stats.std_dev == pytest.approx(1.4142135623730951)
    assert stats.valid_count == 5

    stats = ar.GetStatistics(False, False)
    assert stats is None

    stats = ar.GetStatistics(False, True)
    assert stats is not None
    assert stats.min == 1.0
    assert stats.max == 5.0
    assert stats.mean == 3.0
    assert stats.std_dev == pytest.approx(1.4142135623730951)
    assert stats.valid_count == 5


def test_mem_md_array_statistics_float32():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim0 = rg.CreateDimension("dim0", "unspecified type", "unspecified direction", 2)
    dim1 = rg.CreateDimension("dim1", "unspecified type", "unspecified direction", 3)
    float32dt = gdal.ExtendedDataType.Create(gdal.GDT_Float32)
    ar = rg.CreateMDArray("myarray", [dim0, dim1], float32dt)
    ar.SetNoDataValueDouble(6)
    data = struct.pack("f" * 6, 1, 2, 3, 4, 5, 6)
    ar.Write(data)

    stats = ar.ComputeStatistics(False)
    assert stats.min == 1.0
    assert stats.max == 5.0
    assert stats.mean == 3.0
    assert stats.std_dev == pytest.approx(1.4142135623730951)
    assert stats.valid_count == 5


def test_mem_md_array_copy_autoscale():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim0 = rg.CreateDimension("dim0", "unspecified type", "unspecified direction", 2)
    dim1 = rg.CreateDimension("dim1", "unspecified type", "unspecified direction", 3)
    float32dt = gdal.ExtendedDataType.Create(gdal.GDT_Float32)
    ar = rg.CreateMDArray("myarray", [dim0, dim1], float32dt)
    data = struct.pack("f" * 6, 1.5, 2, 3, 4, 5, 6.5)
    ar.Write(data)
    attr = ar.CreateAttribute(
        "attr_float64", [1], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    attr.Write(1.25)

    attr = ar.CreateAttribute(
        "valid_min", [1], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    attr.Write(1.25)

    out_ds = drv.CreateCopy("", ds, options=["ARRAY:AUTOSCALE=YES"])
    out_rg = out_ds.GetRootGroup()
    out_ar = out_rg.OpenMDArray("myarray")
    assert out_ar.GetAttribute("attr_float64") is not None
    assert out_ar.GetAttribute("valid_min") is None
    assert out_ar.GetDataType() == gdal.ExtendedDataType.Create(gdal.GDT_UInt16)
    assert out_ar.GetOffset() == 1.5
    assert out_ar.GetScale() == (6.5 - 1.5) / 65535.0
    assert struct.unpack("H" * 6, out_ar.Read()) == (
        0,
        6554,
        19661,
        32768,
        45875,
        65535,
    )
    assert struct.unpack("d" * 6, out_ar.GetUnscaled().Read()) == pytest.approx(
        (1.5, 2, 3, 4, 5, 6.5), abs=out_ar.GetScale() / 2
    )


def test_mem_md_array_copy_autoscale_with_explicit_data_type_and_nodata():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim0 = rg.CreateDimension("dim0", "unspecified type", "unspecified direction", 2)
    dim1 = rg.CreateDimension("dim1", "unspecified type", "unspecified direction", 3)
    float32dt = gdal.ExtendedDataType.Create(gdal.GDT_Float32)
    ar = rg.CreateMDArray("myarray", [dim0, dim1], float32dt)
    ar.SetNoDataValueDouble(5)
    data = struct.pack("f" * 6, 1.5, 2, 3, 4, 6.5, 5)
    ar.Write(data)

    out_ds = drv.CreateCopy(
        "", ds, options=["ARRAY:AUTOSCALE=YES", "ARRAY:AUTOSCALE_DATA_TYPE=Int16"]
    )
    out_rg = out_ds.GetRootGroup()
    out_ar = out_rg.OpenMDArray("myarray")
    assert out_ar.GetDataType() == gdal.ExtendedDataType.Create(gdal.GDT_Int16)
    assert out_ar.GetScale() == (6.5 - 1.5) / (65535.0 - 1)
    assert out_ar.GetOffset() == 1.5 - (-32768) * out_ar.GetScale()
    assert out_ar.GetNoDataValueAsDouble() == 32767.0
    assert struct.unpack("h" * 6, out_ar.Read()) == (
        -32768,
        -26215,
        -13108,
        -1,
        32766,
        32767,
    )
    unscaled = struct.unpack("d" * 6, out_ar.GetUnscaled().Read())
    assert unscaled[0:5] == pytest.approx(
        (1.5, 2, 3, 4, 6.5), abs=out_ar.GetScale() / 2
    )
    assert math.isnan(unscaled[5])


def test_mem_md_array_nodata_int64():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim0 = rg.CreateDimension("dim0", "unspecified type", "unspecified direction", 2)
    dim1 = rg.CreateDimension("dim1", "unspecified type", "unspecified direction", 2)
    myarray = rg.CreateMDArray(
        "myarray", [dim0, dim1], gdal.ExtendedDataType.Create(gdal.GDT_Int64)
    )
    assert myarray
    val = -(1 << 63)
    assert myarray.SetNoDataValue(val) == gdal.CE_None
    assert myarray.GetNoDataValue() == val

    # Test MultiDim -> GDALDataset bridge
    ds = myarray.AsClassicDataset(0, 1)
    assert ds.GetRasterBand(1).GetNoDataValue() == val


def test_mem_md_array_nodata_uint64():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim0 = rg.CreateDimension("dim0", "unspecified type", "unspecified direction", 2)
    dim1 = rg.CreateDimension("dim1", "unspecified type", "unspecified direction", 2)
    myarray = rg.CreateMDArray(
        "myarray", [dim0, dim1], gdal.ExtendedDataType.Create(gdal.GDT_UInt64)
    )
    assert myarray
    val = (1 << 64) - 1
    assert myarray.SetNoDataValue(val) == gdal.CE_None
    assert myarray.GetNoDataValue() == val

    # Test MultiDim -> GDALDataset bridge
    ds = myarray.AsClassicDataset(0, 1)
    assert ds.GetRasterBand(1).GetNoDataValue() == val


def test_mem_md_getcoordinatevariables():

    drv = gdal.GetDriverByName("MEM")
    mem_ds = drv.CreateMultiDimensional("myds")
    rg = mem_ds.GetRootGroup()

    dimOther = rg.CreateDimension("other", None, None, 4)
    other = rg.CreateMDArray(
        "other", [dimOther], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    assert other

    dimNode = rg.CreateDimension("node", None, None, 4)
    varX = rg.CreateMDArray(
        "varX", [dimNode], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    assert varX
    varY = rg.CreateMDArray(
        "varY", [dimNode], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    assert varY

    ar = rg.CreateMDArray(
        "ar", [dimOther, dimNode], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    coordinates = ar.CreateAttribute(
        "coordinates", [], gdal.ExtendedDataType.CreateString()
    )
    assert coordinates.WriteString("other varY varX") == 0

    coordVars = ar.GetCoordinateVariables()
    assert len(coordVars) == 3

    assert coordinates.WriteString("other non_existing") == 0
    with gdal.quiet_errors():
        coordVars = ar.GetCoordinateVariables()
        assert len(coordVars) == 1


@pytest.mark.parametrize("new_size", [[], [1, 2], [0]])
def test_mem_md_resize_dim_wrong_new_size(new_size):

    drv = gdal.GetDriverByName("MEM")
    mem_ds = drv.CreateMultiDimensional("myds")
    rg = mem_ds.GetRootGroup()
    xDim = rg.CreateDimension("x", None, None, 4)
    v = rg.CreateMDArray("v", [xDim], gdal.ExtendedDataType.Create(gdal.GDT_Float64))
    v.Write(struct.pack("d" * 4, 1, 2, 3, 4))
    with gdal.quiet_errors():
        assert v.Resize(new_size) == gdal.CE_Failure


def test_mem_md_resize_dim_wrong_too_big_allocation_before_malloc():

    drv = gdal.GetDriverByName("MEM")
    mem_ds = drv.CreateMultiDimensional("myds")
    rg = mem_ds.GetRootGroup()
    xDim = rg.CreateDimension("x", None, None, 4)
    v = rg.CreateMDArray("v", [xDim], gdal.ExtendedDataType.Create(gdal.GDT_Float64))
    v.Write(struct.pack("d" * 4, 1, 2, 3, 4))
    with gdal.quiet_errors():
        assert v.Resize([1 << 63]) == gdal.CE_Failure


@pytest.mark.skipif(sys.maxsize < 2**32, reason="only on 64bit")
def test_mem_md_resize_dim_wrong_too_big_allocation_at_malloc():

    drv = gdal.GetDriverByName("MEM")
    mem_ds = drv.CreateMultiDimensional("myds")
    rg = mem_ds.GetRootGroup()
    xDim = rg.CreateDimension("x", None, None, 4)
    v = rg.CreateMDArray("v", [xDim], gdal.ExtendedDataType.Create(gdal.GDT_Float64))
    v.Write(struct.pack("d" * 4, 1, 2, 3, 4))
    sizeof_double = 8
    with gdal.quiet_errors():
        assert v.Resize([sys.maxsize // sizeof_double]) == gdal.CE_Failure


@pytest.mark.parametrize("new_size,new_values", [(5, (1, 2, 3, 4, 0)), (3, (1, 2, 3))])
def test_mem_md_resize_dim(new_size, new_values):

    drv = gdal.GetDriverByName("MEM")
    mem_ds = drv.CreateMultiDimensional("myds")
    rg = mem_ds.GetRootGroup()
    xDim = rg.CreateDimension("x", None, None, 4)
    v = rg.CreateMDArray("v", [xDim], gdal.ExtendedDataType.Create(gdal.GDT_Float64))
    v.Write(struct.pack("d" * 4, 1, 2, 3, 4))
    assert v.Resize([new_size]) == gdal.CE_None
    assert xDim.GetSize() == new_size
    assert struct.unpack("d" * new_size, v.Read()) == new_values


@pytest.mark.parametrize(
    "new_size,new_values",
    [
        (5, (1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 0.0, 0.0)),
        (3, (1.0, 2.0, 3.0, 4.0, 5.0, 6.0)),
    ],
)
def test_mem_md_resize_first_dim(new_size, new_values):

    drv = gdal.GetDriverByName("MEM")
    mem_ds = drv.CreateMultiDimensional("myds")
    rg = mem_ds.GetRootGroup()
    yDim = rg.CreateDimension("y", None, None, 4)
    xDim = rg.CreateDimension("x", None, None, 2)
    v = rg.CreateMDArray(
        "v", [yDim, xDim], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    v.Write(struct.pack("d" * 8, 1, 2, 3, 4, 5, 6, 7, 8))
    assert v.Resize([new_size, xDim.GetSize()]) == gdal.CE_None
    assert yDim.GetSize() == new_size
    assert xDim.GetSize() == 2
    assert struct.unpack("d" * new_size * 2, v.Read()) == new_values


@pytest.mark.parametrize(
    "new_size,new_values",
    [
        (5, (1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 0.0, 0.0)),
        (3, (1.0, 2.0, 3.0, 4.0, 5.0, 6.0)),
    ],
)
def test_mem_md_resize_first_dim_and_other_array(new_size, new_values):

    drv = gdal.GetDriverByName("MEM")
    mem_ds = drv.CreateMultiDimensional("myds")
    rg = mem_ds.GetRootGroup()
    yDim = rg.CreateDimension("y", None, None, 4)
    xDim = rg.CreateDimension("x", None, None, 2)
    v = rg.CreateMDArray(
        "v", [yDim, xDim], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    v.Write(struct.pack("d" * 8, 1, 2, 3, 4, 5, 6, 7, 8))
    v2 = rg.CreateMDArray(
        "v2", [yDim, xDim], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    v.Write(struct.pack("d" * 8, 1, 2, 3, 4, 5, 6, 7, 8))
    v2.Write(struct.pack("d" * 8, 1, 2, 3, 4, 5, 6, 7, 8))
    assert v.Resize([new_size, xDim.GetSize()]) == gdal.CE_None
    assert yDim.GetSize() == new_size
    assert xDim.GetSize() == 2
    assert struct.unpack("d" * new_size * 2, v.Read()) == new_values
    assert struct.unpack("d" * new_size * 2, v2.Read()) == new_values


@pytest.mark.parametrize(
    "new_size,new_values",
    [
        (5, (1.0, 2.0, 3.0, 4.0, 5.0, 6.0, 7.0, 8.0, 0.0, 0.0)),
        (3, (1.0, 2.0, 3.0, 4.0, 5.0, 6.0)),
    ],
)
def test_mem_md_resize_arbitrary_dim(new_size, new_values):

    drv = gdal.GetDriverByName("MEM")
    mem_ds = drv.CreateMultiDimensional("myds")
    rg = mem_ds.GetRootGroup()
    yDim = rg.CreateDimension("y", None, None, 4)
    xDim = rg.CreateDimension("x", None, None, 2)
    otherDim = rg.CreateDimension("otherDim", None, None, 1)
    v = rg.CreateMDArray(
        "v", [otherDim, yDim, xDim], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    v.Write(struct.pack("d" * 8, 1, 2, 3, 4, 5, 6, 7, 8))
    v2 = rg.CreateMDArray(
        "v2", [otherDim, yDim, xDim], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    v.Write(struct.pack("d" * 8, 1, 2, 3, 4, 5, 6, 7, 8))
    v2.Write(struct.pack("d" * 8, 1, 2, 3, 4, 5, 6, 7, 8))

    with gdal.quiet_errors():
        assert v.Resize([1, 1 << 63, xDim.GetSize()]) == gdal.CE_Failure

    assert v.Resize([1, new_size, xDim.GetSize()]) == gdal.CE_None
    assert otherDim.GetSize() == 1
    assert yDim.GetSize() == new_size
    assert xDim.GetSize() == 2
    assert struct.unpack("d" * new_size * 2, v.Read()) == new_values
    assert struct.unpack("d" * new_size * 2, v2.Read()) == new_values


def test_mem_md_resize_two_dims_at_once():

    drv = gdal.GetDriverByName("MEM")
    mem_ds = drv.CreateMultiDimensional("myds")
    rg = mem_ds.GetRootGroup()
    yDim = rg.CreateDimension("y", None, None, 4)
    xDim = rg.CreateDimension("x", None, None, 2)
    v = rg.CreateMDArray(
        "v", [yDim, xDim], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    v.Write(struct.pack("d" * 8, 1, 2, 3, 4, 5, 6, 7, 8))
    assert v.Resize([3, 4]) == gdal.CE_None
    assert struct.unpack("d" * 3 * 4, v.Read()) == (
        1.0,
        2.0,
        0.0,
        0.0,
        3.0,
        4.0,
        0.0,
        0.0,
        5.0,
        6.0,
        0.0,
        0.0,
    )


def test_mem_md_resize_dim_referenced_twice():

    drv = gdal.GetDriverByName("MEM")
    mem_ds = drv.CreateMultiDimensional("myds")
    rg = mem_ds.GetRootGroup()
    xDim = rg.CreateDimension("x", None, None, 2)
    v = rg.CreateMDArray(
        "v", [xDim, xDim], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    v.Write(struct.pack("d" * 4, 1, 2, 3, 4))
    assert v.Resize([3, 3]) == gdal.CE_None
    assert struct.unpack("d" * 3 * 3, v.Read()) == (
        1.0,
        2.0,
        0.0,
        3.0,
        4.0,
        0.0,
        0.0,
        0.0,
        0.0,
    )


def test_mem_md_resize_dim_referenced_twice_error():

    drv = gdal.GetDriverByName("MEM")
    mem_ds = drv.CreateMultiDimensional("myds")
    rg = mem_ds.GetRootGroup()
    xDim = rg.CreateDimension("x", None, None, 2)
    v = rg.CreateMDArray(
        "v", [xDim, xDim], gdal.ExtendedDataType.Create(gdal.GDT_Float64)
    )
    with gdal.quiet_errors():
        assert v.Resize([2, 3]) == gdal.CE_Failure
        assert v.Resize([3, 2]) == gdal.CE_Failure


@gdaltest.enable_exceptions()
def test_mem_md_rename_group():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    group = rg.CreateGroup("group")
    group_attr = group.CreateAttribute(
        "group_attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    rg.CreateGroup("other_group")
    dim = group.CreateDimension("dim0", "unspecified type", "unspecified direction", 2)
    ar = group.CreateMDArray("ar", [dim], gdal.ExtendedDataType.Create(gdal.GDT_Byte))
    ar_attr = ar.CreateAttribute(
        "ar_attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )

    subgroup = group.CreateGroup("subgroup")
    subgroup.CreateGroup("subsubgroup")
    subgroup_attr = subgroup.CreateAttribute(
        "subgroup_attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    subgroup_ar = subgroup.CreateMDArray(
        "subgroup_ar", [dim], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    subgroup_ar_attr = subgroup_ar.CreateAttribute(
        "subgroup_ar_attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )

    # Cannot rename root group
    with pytest.raises(Exception):
        rg.Rename("foo")

    # Empty name
    with pytest.raises(Exception):
        group.Rename("")

    # Existing name
    with pytest.raises(Exception):
        group.Rename("other_group")

    # Rename group and test effects
    group.Rename("group_renamed")
    assert group.GetName() == "group_renamed"
    assert group.GetFullName() == "/group_renamed"

    assert set(rg.GetGroupNames()) == {"group_renamed", "other_group"}

    assert dim.GetName() == "dim0"
    assert dim.GetFullName() == "/group_renamed/dim0"

    assert group_attr.GetName() == "group_attr"
    assert group_attr.GetFullName() == "/group_renamed/group_attr"

    assert ar.GetName() == "ar"
    assert ar.GetFullName() == "/group_renamed/ar"

    assert ar_attr.GetName() == "ar_attr"
    assert ar_attr.GetFullName() == "/group_renamed/ar/ar_attr"

    assert subgroup.GetName() == "subgroup"
    assert subgroup.GetFullName() == "/group_renamed/subgroup"

    assert subgroup_attr.GetName() == "subgroup_attr"
    assert subgroup_attr.GetFullName() == "/group_renamed/subgroup/subgroup_attr"

    assert subgroup_ar.GetName() == "subgroup_ar"
    assert subgroup_ar.GetFullName() == "/group_renamed/subgroup/subgroup_ar"

    assert subgroup_ar_attr.GetName() == "subgroup_ar_attr"
    assert (
        subgroup_ar_attr.GetFullName()
        == "/group_renamed/subgroup/subgroup_ar/subgroup_ar_attr"
    )


@gdaltest.enable_exceptions()
def test_mem_md_rename_dimension():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    dim = rg.CreateDimension("dim", "unspecified type", "unspecified direction", 2)
    rg.CreateDimension("other_dim", "unspecified type", "unspecified direction", 2)
    ar = rg.CreateMDArray("ar", [dim], gdal.ExtendedDataType.Create(gdal.GDT_Byte))

    # Empty name
    with pytest.raises(Exception):
        dim.Rename("")

    # Existing name
    with pytest.raises(Exception):
        dim.Rename("other_dim")

    # Rename dimension and test effects
    dim.Rename("dim_renamed")
    assert dim.GetName() == "dim_renamed"
    assert dim.GetFullName() == "/dim_renamed"

    assert set(x.GetFullName() for x in rg.GetDimensions()) == {
        "/dim_renamed",
        "/other_dim",
    }

    assert [x.GetFullName() for x in ar.GetDimensions()] == ["/dim_renamed"]


@gdaltest.enable_exceptions()
def test_mem_md_rename_attribute():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    subg = rg.CreateGroup("group")
    subg_attr = subg.CreateAttribute(
        "subg_attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    subg.CreateAttribute(
        "subg_other_attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    ar = subg.CreateMDArray("ar", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte))
    attr = ar.CreateAttribute("attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte))
    ar.CreateAttribute("other_attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte))

    # Empty name
    with pytest.raises(Exception):
        attr.Rename("")

    # Existing name
    with pytest.raises(Exception):
        attr.Rename("other_attr")

    # Rename array attribute and test effects
    attr.Rename("attr_renamed")
    assert attr.GetName() == "attr_renamed"
    assert attr.GetFullName() == "/group/ar/attr_renamed"
    assert set(x.GetName() for x in ar.GetAttributes()) == {
        "attr_renamed",
        "other_attr",
    }

    # Existing name
    with pytest.raises(Exception):
        subg_attr.Rename("subg_other_attr")

    # Rename group attribute and test effects
    subg_attr.Rename("subg_attr_renamed")
    assert subg_attr.GetName() == "subg_attr_renamed"
    assert subg_attr.GetFullName() == "/group/_GLOBAL_/subg_attr_renamed"
    assert set(x.GetName() for x in subg.GetAttributes()) == {
        "subg_attr_renamed",
        "subg_other_attr",
    }


@gdaltest.enable_exceptions()
def test_mem_md_rename_array():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()
    subg = rg.CreateGroup("group")
    ar = subg.CreateMDArray("ar", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte))
    subg.CreateMDArray("other_array", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte))
    attr = ar.CreateAttribute("attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte))

    # Empty name
    with pytest.raises(Exception):
        ar.Rename("")

    # Existing name
    with pytest.raises(Exception):
        ar.Rename("other_array")

    # Rename array and test effects
    ar.Rename("ar_renamed")
    assert ar.GetName() == "ar_renamed"
    assert ar.GetFullName() == "/group/ar_renamed"

    assert attr.GetName() == "attr"
    assert attr.GetFullName() == "/group/ar_renamed/attr"

    assert set(subg.GetMDArrayNames()) == {"ar_renamed", "other_array"}


@gdaltest.enable_exceptions()
def test_mem_md_delete_group():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    group = rg.CreateGroup("group")
    group_attr = group.CreateAttribute(
        "ar_attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    ar = group.CreateMDArray("ar", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte))
    ar_attr = ar.CreateAttribute(
        "ar_attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )
    subgroup = group.CreateGroup("subgroup")

    with pytest.raises(Exception, match="is not a sub-group"):
        rg.DeleteGroup("non_existing")

    # Delete group and test effects
    rg.DeleteGroup("group")

    assert len(rg.GetGroupNames()) == 0

    with pytest.raises(Exception, match="has been deleted"):
        group.GetMDArrayNames()

    with pytest.raises(Exception, match="has been deleted"):
        group_attr.Rename("foo")

    with pytest.raises(Exception, match="has been deleted"):
        subgroup.GetMDArrayNames()

    with pytest.raises(Exception, match="has been deleted"):
        ar.GetAttributes()

    with pytest.raises(Exception, match="has been deleted"):
        ar_attr.Rename("foo")


@gdaltest.enable_exceptions()
def test_mem_md_delete_array():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    group = rg.CreateGroup("group")
    ar = group.CreateMDArray("ar", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte))
    ar_attr = ar.CreateAttribute(
        "ar_attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )

    with pytest.raises(Exception, match="is not an array"):
        group.DeleteMDArray("non_existing")

    # Delete array and test effects
    group.DeleteMDArray("ar")

    assert len(rg.GetMDArrayNames()) == 0

    with pytest.raises(Exception, match="has been deleted"):
        ar.GetAttributes()

    with pytest.raises(Exception, match="has been deleted"):
        ar_attr.Rename("foo")


@gdaltest.enable_exceptions()
def test_mem_md_delete_group_attribute():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    group = rg.CreateGroup("group")
    group_attr = group.CreateAttribute(
        "group_attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )

    with pytest.raises(Exception, match="is not an attribute"):
        group.DeleteAttribute("non_existing")

    # Delete attribute and test effects
    group.DeleteAttribute("group_attr")

    assert len(group.GetAttributes()) == 0

    with pytest.raises(Exception, match="has been deleted"):
        group_attr.Rename("foo")


@gdaltest.enable_exceptions()
def test_mem_md_delete_array_attribute():

    drv = gdal.GetDriverByName("MEM")
    ds = drv.CreateMultiDimensional("myds")
    rg = ds.GetRootGroup()

    group = rg.CreateGroup("group")
    ar = group.CreateMDArray("ar", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte))
    ar_attr = ar.CreateAttribute(
        "ar_attr", [], gdal.ExtendedDataType.Create(gdal.GDT_Byte)
    )

    with pytest.raises(Exception, match="is not an attribute"):
        ar.DeleteAttribute("non_existing")

    # Delete attribute and test effects
    ar.DeleteAttribute("ar_attr")

    assert len(ar.GetAttributes()) == 0

    with pytest.raises(Exception, match="has been deleted"):
        ar_attr.Rename("foo")


def XX_test_all_forever():
    while True:
        test_mem_md_basic()
        test_mem_md_subgroup()
        test_mem_md_array_unnamed_array()
        test_mem_md_array_duplicated_array_name()
        test_mem_md_array_nodim()
        test_mem_md_array_single_dim()
        test_mem_md_array_string()
        test_mem_md_datatypes()
        test_mem_md_array_compoundtype()
        test_mem_md_array_3_dim()
        test_mem_md_array_4_dim()
        test_mem_md_copy_array()
        test_mem_md_array_read_write_errors()
        test_mem_md_invalid_dims()
        test_mem_md_array_invalid_args()
        test_mem_md_array_too_large()
        test_mem_md_array_too_large_overflow_dim()
        test_mem_md_array_30dim()
        test_mem_md_array_32dim()
        test_mem_md_group_attribute_single_numeric()
        test_mem_md_group_attribute_multiple_numeric()
        test_mem_md_group_attribute_single_string()
        test_mem_md_group_attribute_multiple_string()
        test_mem_md_array_attribute()
        test_mem_md_array_slice()
        test_mem_md_band_as_mdarray()
        test_mem_md_array_as_classic_dataset()
        test_mem_md_array_transpose()
        test_mem_md_array_single_dim_non_contiguous_copy()