Creation of mesh file with dataset for UGRID

[jmkerloch]

I created a repository to test mesh file creation with mdal-python : https://gitlab.com/jmkerloch/test-mdal-python.

I'm trying to create mesh file with all available drivers from my MDAL build.

I managed to create some mesh file but I don't understand why there is no data when using UGRID driver.

from mdal import PyMesh, drivers, MDAL_DataLocation
import numpy as np

def get_ugrid_driver():
    for test_driver in drivers():
        if test_driver.name == "Ugrid":
            return test_driver
    return None

driver = get_ugrid_driver()
test = PyMesh(driver)

# Create vertices and faces
test.vertices=np.array([(0. , 0., 0.), (1. , 0., 0.), (2. , 0., 0.), (0.5, 1., 0.), (1.5, 1., 0.), (1. , 2., 0.)], dtype=[('X', '<f8'), ('Y', '<f8'), ('Z', '<f8')])
test.faces=np.array([(3, 0, 1, 3), (3, 1, 4, 3), (3, 1, 2, 4), (3, 3, 4, 5)], dtype=[('Vertices', '<u4'), ('V0', '<u4'), ('V1', '<u4'), ('V2', '<u4')])

# Add group for data on faces
group= test.add_group("test_data_on_face", location=MDAL_DataLocation.DataOnFaces)

# Define data  
faces_val = [(i,) for i in range(test.faces.shape[0])]
data = np.array(faces_val, dtype=[('U', '<f8')])

# Add to group
group.add_data(data)    

test.save(f"created_ugrid_mesh.nc")

The mesh is displayed in QGIS but there is no dataset test_data_on_face

@runette could you provide some guidance for ugrid mesh with dataset creation ?

[jmkerloch]

@runette It seems that we need to indicate the uri of the group so we can save it. In my code there was also missing the save of the group.

I managed to create a Ugrid file with time support. Here is the python code.

from mdal import PyMesh, drivers, MDAL_DataLocation
import numpy as np

def get_ugrid_driver():
    for test_driver in drivers():
        if test_driver.name == "Ugrid":
            return test_driver
    return None

driver = get_ugrid_driver()
test = PyMesh(driver)

# Create vertices and faces
test.vertices=np.array([(0. , 0., 0.), (1. , 0., 0.), (2. , 0., 0.), (0.5, 1., 0.), (1.5, 1., 0.), (1. , 2., 0.)], dtype=[('X', '<f8'), ('Y', '<f8'), ('Z', '<f8')])
test.faces=np.array([(3, 0, 1, 3), (3, 1, 4, 3), (3, 1, 2, 4), (3, 3, 4, 5)], dtype=[('Vertices', '<u4'), ('V0', '<u4'), ('V1', '<u4'), ('V2', '<u4')])

# Add group for data on faces
group= test.add_group("test_data_on_face", location=MDAL_DataLocation.DataOnFaces, uri="Ugrid:created_ugrid_mesh.nc")

# Define data  
faces_val = [(i,) for i in range(test.faces.shape[0])]
data = np.array(faces_val, dtype=[('U', '<f8')])

reversed_faces_val = [(i,) for i in reversed(range(test.faces.shape[0]))]
reversed_data = np.array(reversed_faces_val, dtype=[('U', '<f8')])

# Add to group
nb_time = 24
for time in range(nb_time):
    if time%2:
        group.add_data(data, float(time))
    else:
        group.add_data(reversed_data, float(time))

group.save()

test.save(f"created_ugrid_mesh.nc") 

For me it seems that there are no issue with mdal-python. I'm closing this issue.

[runette]

Sorry for the delay - I have been traveling

I do really need to create better documentation - I know.

The control flow is based on the MDAL C API to keep the bindings fairly simple.

However, I would use a slightly different flow. The reason you are having to give a group URI is that the mesh does not yet have a uri - this is a requirement of the C API

I would do something like ...

ds = Datasource("Ugrid:created_ugrid_mesh.nc")

With ds.add_mesh() as test: .... Create the vertices and faces etc test.add_group( ... ) # should not need a uri test.save()

Etc - should work but I am currently on my phone so cannot test it and have probably made some typos :)

On Tue, 30 May 2023, 10:36 jmkerloch, @.***> wrote:

@runette https://github.com/runette It seems that we need to indicate the uri of the group so we can save it. In my code there was also missing the save of the group.

I managed to create a Ugrid file with time support. Here is the python code.

from mdal import PyMesh, drivers, MDAL_DataLocationimport numpy as np def get_ugrid_driver(): for test_driver in drivers(): if test_driver.name == "Ugrid": return test_driver return None driver = get_ugrid_driver()test = PyMesh(driver)

Create vertices and facestest.vertices=np.array([(0. , 0., 0.), (1. , 0., 0.), (2. , 0., 0.), (0.5, 1., 0.), (1.5, 1., 0.), (1. , 2., 0.)], dtype=[('X', '<f8'), ('Y', '<f8'), ('Z', '<f8')])test.faces=np.array([(3, 0, 1, 3), (3, 1, 4, 3), (3, 1, 2, 4), (3, 3, 4, 5)], dtype=[('Vertices', '<u4'), ('V0', '<u4'), ('V1', '<u4'), ('V2', '<u4')])

Add group for data on facesgroup= test.add_group("test_data_on_face", location=MDAL_DataLocation.DataOnFaces, uri="Ugrid:created_ugrid_mesh.nc")

Define data faces_val = [(i,) for i in range(test.faces.shape[0])]data = np.array(faces_val, dtype=[('U', '<f8')])

reversed_faces_val = [(i,) for i in reversed(range(test.faces.shape[0]))]reversed_data = np.array(reversed_faces_val, dtype=[('U', '<f8')])

Add to groupnb_time = 24for time in range(nb_time):

if time%2:
    group.add_data(data, float(time))
else:
    group.add_data(reversed_data, float(time))

group.save() test.save(f"created_ugrid_mesh.nc")

For me it seems that there are no issue with mdal-python. I'm closing this issue.

— Reply to this email directly, view it on GitHub https://github.com/ViRGIS-Team/mdal-python/issues/8#issuecomment-1568006848, or unsubscribe https://github.com/notifications/unsubscribe-auth/AARC2MZYN5P6YOCJKZXHA3LXIWWSPANCNFSM6AAAAAAYS5HJH4 . You are receiving this because you were mentioned.Message ID: @.***>

[jmkerloch]

I tried to add the uri to PyMesh but I must use a Datasource.

I was trying wihout a datasource because I will only have one mesh on my file.

[runette]

Both ways work.

I just think the datasource way is more intuitive and it also automatically closes the PyMesh object - which is a wrapper around the C object - cleanly.

On Wed, 31 May 2023, 12:17 jmkerloch, @.***> wrote:

I tried to add the uri to PyMesh but I must use a Datasource.

I was trying wihout a datasource because I will only have one mesh on my file.

— Reply to this email directly, view it on GitHub https://github.com/ViRGIS-Team/mdal-python/issues/8#issuecomment-1569912935, or unsubscribe https://github.com/notifications/unsubscribe-auth/AARC2M5KXGZZ6VS2XYVFGBDXI4LEZANCNFSM6AAAAAAYS5HJH4 . You are receiving this because you were mentioned.Message ID: @.***>

[runette]

This worked for me

from mdal import Datasource, MDAL_DataLocation import numpy as np

ds = Datasource("created_ugrid_mesh.nc", "Ugrid")

with ds.add_mesh() as test: test.vertices=np.array([(0. , 0., 0.), (1. , 0., 0.), (2. , 0., 0.), (0.5, 1., 0.), (1.5, 1., 0.), (1. , 2., 0.)], dtype=[('X', '<f8'), ('Y', '<f8'), ('Z', '<f8')]) test.faces=np.array([(3, 0, 1, 3), (3, 1, 4, 3), (3, 1, 2, 4), (3, 3, 4, 5)], dtype=[('Vertices', '<u4'), ('V0', '<u4'), ('V1', '<u4'), ('V2', '<u4')])

Add group for data on faces

group= test.add_group("test_data_on_face", location= MDAL_DataLocation.DataOnFaces)

Define data

faces_val = [(i,) for i in range(test.faces.shape[0])] data = np.array(faces_val, dtype=[('U', '<f8')])

reversed_faces_val = [(i,) for i in reversed(range(test.faces.shape[0]))] reversed_data = np.array(reversed_faces_val, dtype=[('U', '<f8')])

Add to group

nb_time = 24 for time in range(nb_time): if time%2: group.add_data(data, float(time)) else: group.add_data(reversed_data, float(time)) group.save() test.save()

On Wed, 31 May 2023 at 12:23, Paul Harwood @.***> wrote:

Both ways work.

I just think the datasource way is more intuitive and it also automatically closes the PyMesh object - which is a wrapper around the C object - cleanly.

On Wed, 31 May 2023, 12:17 jmkerloch, @.***> wrote:

I tried to add the uri to PyMesh but I must use a Datasource.

I was trying wihout a datasource because I will only have one mesh on my file.

— Reply to this email directly, view it on GitHub https://github.com/ViRGIS-Team/mdal-python/issues/8#issuecomment-1569912935, or unsubscribe https://github.com/notifications/unsubscribe-auth/AARC2M5KXGZZ6VS2XYVFGBDXI4LEZANCNFSM6AAAAAAYS5HJH4 . You are receiving this because you were mentioned.Message ID: @.***>

[runette]

That got mangled - should have been

from mdal import Datasource, MDAL_DataLocation
import numpy as np

ds = Datasource("created_ugrid_mesh.nc", "Ugrid")

with ds.add_mesh() as test:
    test.vertices=np.array([(0. , 0., 0.), (1. , 0., 0.), (2. , 0., 0.), (0.5, 1., 0.), (1.5, 1., 0.), (1. , 2., 0.)], dtype=[('X', '<f8'), ('Y', '<f8'), ('Z', '<f8')])
    test.faces=np.array([(3, 0, 1, 3), (3, 1, 4, 3), (3, 1, 2, 4), (3, 3, 4, 5)], dtype=[('Vertices', '<u4'), ('V0', '<u4'), ('V1', '<u4'), ('V2', '<u4')])

# Add group for data on faces
    group= test.add_group("test_data_on_face", location=MDAL_DataLocation.DataOnFaces)

# Define data  
    faces_val = [(i,) for i in range(test.faces.shape[0])]
    data = np.array(faces_val, dtype=[('U', '<f8')])

    reversed_faces_val = [(i,) for i in reversed(range(test.faces.shape[0]))]
    reversed_data = np.array(reversed_faces_val, dtype=[('U', '<f8')])

# Add to group
    nb_time = 24
    for time in range(nb_time):
        if time%2:
            group.add_data(data, float(time))
        else:
            group.add_data(reversed_data, float(time))
    group.save()
    test.save()