Processing .asc Point Clouds/ .stl Mesh files into 3D Integer Numpy Arrays and Indexing-based material's Assigment

Hello dear community, fellow pro developers and Python helpers!

I hope everyone is having a nice day.
I come to you with a big challenge! I am currently working in processing 3D CAD models via .stl mesh files into Python. I am working with an EM simulator that uses .HDF5 files to insert geometry. The purpose is to export multiple meshes from a FreeCAD file that make up the entire compound and convert them into a 3D Integer-indexed materials NumPy array. The Python lib takes as parameters the .stl of the feature part bodies of the 3D CAD and converts them via this .HDF5 file. I attach to you the script of the aforementioned library(convert_mesh_3DNumpy.py). Part bodies of the CAD compound are extracted using the pd_extract_features.py Macro that simply extracts each part body from the compound and shows them individually in the FreeCAD interface.

Further work and efforts are to be put in order to integrate multiple .stl files (ie .stl file List in Python) to convert all meshes.
This CAD model is employed in EM applications. For each part (feature) it is then needed to add electric material properties (relative permittivity, conductivity, permeability, and magnetic loss) to each one individually.

The simulator that gets the Python lib output then uses an insert geometry command that takes the .HDF5 integer-based file and a .txt text file materials. An example of the material’s file is attached.

  1. So my question is which CAD tool WB can use the Points WB, Mesh WB or whatever tool is the most effective to properly convert all the part body (ie, features) that make up a body in Python Numpy grid??

  2. The second problem is how can I use a macro, tools, Python lib, whatever it is that can allow me to assign these electric properties aforementioned from FreeCAD 3D model directly into a .txt text file in the format attached, preferably via the Mesh-file processing?

I would very much appreciate it if you amazing and grateful community can advise me in any way!

Kind regards,

Diego

Attached files:

  1. binary_voxel_lib.py (convert .asc Point Clouds to Binary Voxel geometry files)

  2. convert_mesh_3DNumpy( convert .stl Mesh files to Binary Voxel geometry files)

  3. Materials .txt file format (ie, format that needs to be assigned for each mesh file that makes up the 3D CAD file

  4. Input file .in example of the EM simulator, that takes the .HDF5 geometry file and the .txt Material file into account.

  5. binary_voxel_lib.py (convert .asc Point Clouds to Binary Voxel geometry files)

import sys

import numpy as np
import pandas as pd
import os
from pyntcloud import PyntCloud
import h5py

def convert_points_binaryvoxel(filename_asc, coords_init ,domain, discretization):

"""
Método que convierte el archivo de punto .asc Point Cloud de FreeCAD en binary voxel librería de 0-1.

Args
:param filename_asc: (string) Path - nombre del archivo de puntos exportados de FreeCAD en formato .asc
:param coords_init: (list) Lista de las coordenadas iniciales de cuerpo-body exportado de FreeCAD.
:param domain: (list) Lista de las dimensiones del modelo 3D a simular en gprMax, unidades en m.
:param discretization: (list) Lista de las dimensiones de la resolución espacial del modelo 3D a simular en gprMax.
La resolución debe ser la misma dada por parámetro en el archivo de puntos exportado de FreeCAD.


:return: binary_voxel: (ndarray) Numpy Data Array datos binarios del archivo de puntos. Asigna 0 para coordenadas con espacio libre y 1 para coordenadas con material.
:return: binary_coords: (list) Lista de coordenadas donde hay material dentro del binary voxel libreria
"""

# Abre el archivo
if not os.path.exists(filename_asc):
    sys.exit("Path Error: El path del archivo no existe!!")
else:
    pass



# LOS SEGMENTOS DE DIVISION CAMBIAN PARA CADA PARTE

#segmentos de división
x_seg = int(float(domain[0]/discretization[0]))
y_seg = int(float(domain[1]/discretization[1]))
z_seg = int(float(domain[2]/discretization[2]))

x_delta = coords_init[0]*1000
y_delta = coords_init[1]*1000
z_delta = coords_init[2]*1000

# Lista de coordenadas de material de la binary voxel libreria
binary_coords = []

# Lista de puntos desplazados del archivo de puntos exportados
delta_points = []

# 3D Numpy Array datos con dimensiones de los segmentos del modelo 3D de gprMax
datos = np.zeros((int(domain[0]/discretization[0]), int(domain[1]/discretization[1]), int(domain[2]/discretization[2])), dtype=np.int16)

# carga los puntos Point Cloud del archivo .asc
point_cloud= np.loadtxt(filename_asc)
# le da shape al Point Cloud
point_cloud.shape

# for point in point_cloud:
#     actual_point = point
#
#     x_new = actual_point[0] + x_delta
#     y_new = actual_point[1] + y_delta
#     z_new = actual_point[2] + z_delta
#
#     new_point=[x_new,y_new,z_new]
#
#     delta_points.append(new_point)


# crea el Data Frame con los puntos del Point Cloud y las columnas x,y,z
df=pd.DataFrame(data=point_cloud, columns=['x','y','z'])
# crea el Cloud de nube de puntos PyntCloud
cloud=PyntCloud(df)
# grafica el plot de los puntos Cloud Plot
cloud.plot(initial_point_size=0.1, opacity=0.5)

# agrega el voxel id con la estructura de voxelgrid y los segmentos de división en x,y,z
voxel_id=cloud.add_structure("voxelgrid", n_x=x_seg, n_y=y_seg, n_z=z_seg)

voxel=cloud.structures[voxel_id]
voxel.plot(d=3, mode="density", cmap="hsv")

# Crear (Binary) Numpy Array
binary_voxel=voxel.get_feature_vector(mode='binary')
binary_voxel.shape


for ii in range(int(float(datos.shape[0]))):
    for jj in range(int(float(datos.shape[1]))):
        for kk in range(int(float(datos.shape[2]))):
            actual_coords=[ii,jj,kk]
            if binary_voxel[ii,jj,kk]==1:
                datos[ii,jj,kk] = 1
                binary_coords.append(actual_coords)
            else:
                pass


# Creación de archivo de geometría (HDF5)
hdf5file = os.path.splitext(filename_asc)[0] + '.h5'
print("Archivo HDF5 de geometria creado.")

with h5py.File(hdf5file, 'w') as fout:

    fout.attrs['dx_dy_dz'] = discretization
    fout.create_dataset('data', data=datos)

print("El archivo de geometría Binario tipo Voxel es:")
print(datos)
print(type(binary_voxel))

print("Las coordenadas de material de la libreria binary voxel son:")
print(binary_coords)

return datos, binary_coords
  1. convert_mesh_3DNumpy( convert .stl Mesh files to Binary Voxel geometry files)

import sys
from builtins import dict

import numpy as np
import pandas as pd
import os
from pyntcloud import PyntCloud
import h5py
import trimesh
from stl import mesh
from mpl_toolkits import mplot3d
from matplotlib import pyplot

def convert_mesh_3DNumpy(filename_stl,coords_init,domain,discretization):

"""

Args:
:param filename_stl: (string) Path nombre del archivo .stl mesh tipo malla exportado del modelo 3D de FreeCAD.
:param coords_init: (list) Lista de coordenadas iniciales del archivo 3D  CAD exportado.
:param domain: (list) Lista de las dimensiones del dominio a simular del modelo 3D en gprMax en m.
:param discretization: (list) Lista de las dimensiones de la resolución espacial del modelo 3D a simular en gprMax.
El parámetro de discretización espacial debe ser el mismo que la resolución de la malla exportada de FreeCAD.

:return:
"""

# Abre el archivo
if not os.path.exists(filename_stl):
    sys.exit("Path Error: El archivo de binary voxel libreria no existe!!")
else:
    pass

# segmentos de división
x_seg = int(float(domain[0] / discretization[0]))
y_seg = int(float(domain[1] / discretization[1]))
z_seg = int(float(domain[2] / discretization[2]))

# delta de posiciones iniciales en mm
x_delta = coords_init[0]*1000
y_delta = coords_init[1]*1000
z_delta = coords_init[2]*1000


# 3D Numpy Array datos con dimensiones de los segmentos del modelo 3D de gprMax
datos = np.zeros((int(domain[0]/discretization[0]), int(domain[1]/discretization[1]), int(domain[2]/discretization[2])), dtype=np.int16)

# Point Cloud de puntos extraidos de mesh files
point_cloud = []

# Delta Point de puntos desplazados en las coordenadas iniciales del body

delta_points = []

# carga la malla trimesh del archivo .stl
mesh_tri = trimesh.load_mesh(filename_stl)
assert (mesh_tri.is_watertight)
mesh_tri.show()

# genera voxel grid de la malla mesh con tamaño de voxel de 1mm
vol = mesh_tri.voxelized(pitch=1)
# extrae los puntos de la malla voxelizada
point_cloud = vol.points
point_cloud.shape

# for point in point_cloud:
#     actual_point = point
#
#     x_new = actual_point[0] + x_delta
#     y_new = actual_point[1] + y_delta
#     z_new = actual_point[2] + z_delta
#
#     new_point=[x_new,y_new,z_new]
#
#     delta_points.append(new_point)


# crea el data frame de los puntos extraidos
df = pd.DataFrame(data=point_cloud, columns=['x', 'y', 'z'])
# crea clase Point Cloud nube de puntos con los puntos de malla
cloud = PyntCloud(df)
cloud.plot(initial_point_size=0.1)

voxel_id=cloud.add_structure("voxelgrid", n_x=x_seg, n_y=y_seg, n_z=z_seg)
voxel=cloud.structures[voxel_id]
voxel.plot(d=3, mode="density", cmap="hsv")

# genera 3D Numpy Binary Voxel array a partir de Point Cloud extraido de la malla
binary_voxel=voxel.get_feature_vector(mode='binary')
binary_voxel.shape



# crear mesh malla a partir de archivo .stl exportado de FreeCAD

#mesh_init = mesh.Mesh.from_file(filename_stl)
#print(type(mesh_init))
#mesh_1 = o3d.io.read_triangle_mesh(filename_stl)

#df_init = pd.DataFrame(data=mesh_1, dtype=np.meshgrid, columns=['x', 'y', 'z'])
#cloud_2 = PyntCloud(points=df, mesh=mesh_1)


# recorre la matriz de datos de salida y agrega los datos de 3D Binary Voxel lib
for ii in range(datos.shape[0]):
    for jj in range(datos.shape[1]):
        for kk in range(datos.shape[2]):

            actual_vox = binary_voxel[ii,jj,kk]

            if actual_vox == 1:
                datos[ii,jj,kk] = 1
            else:
                datos[ii, jj, kk] = 0


# Creación de archivo de geometría (HDF5)
hdf5file = os.path.splitext(filename_stl)[0] + '.h5'
print("Archivo HDF5 de geometria creado.")

with h5py.File(hdf5file, 'w') as fout:

    fout.attrs['dx_dy_dz'] = discretization
    fout.create_dataset('data', data=binary_voxel, dtype=np.int16)

#print("Los puntos de malla son:")
#print(delta_points)

print("Los datos del archivo HDF5 de geometria son:")
print(datos)
return datos
  1. Materials .txt file format (ie, format that needs to be assigned for each mesh file that makes up the 3D CAD file

Vivaldi Antenna Material properties

Non-dispersive properties from Dielectric Properties of Body Tissues: HTML clients

#material: 1 0 1 0 free_space
#material: 1 59600000 1 0 Copper1
#material: 2.2 0 1 1 duroid_5880
#material: 1 59600000 1 0 Copper2

  1. Input file .in example of the EM simulator, that takes the .HDF5 geometry file and the .txt Material file into account.

#title: vivaldi_001

#domain: 0.1 0.1 0.05
#dx_dy_dz: 0.0015 0.0015 0.0015
#time_window: 7e-10

LA FRECUENCIA CENTRAL (DE EXCITACIÓN) ES LA FRECUENCIA DE LA FORMA DE ONDA ESPERADA DE LA FUENTE SOURCE

LA FRECUENCIA MAXIMA EN MAXIMA PERMITIVIDAD OCURRE EN EL LAMBDA MÁS PEQUEÑO ( PARA EL MATERIAL DE PERMITIVIDAD MAYOR)

LA FRECUENCIA MÁXIMA ES ENTONCES 3-4 FC VECES LA FRECUENCIA CENTRAL

FUNCIONA CON dx_dy_dz: 0.0015 0.0015 0.0015 (d=0.5mm.asc)

##waveform: ricker 1 10e9 my_ricker
##hertzian_dipole: z 0.100 0.170 0 my_ricker
##rx: 0.140 0.170 0

#python:
from user_libs.construccion.binary_voxel_lib_resp import convert_points_binaryvoxel
convert_points_binaryvoxel(“C:/Free_CAD/WORKSPACE/Export_Files/vivaldi_antenna_Final/d=0.5mm/antenna_final_d=0.5.asc”,[0,0,0.00708],[0.1,0.1,0.05],[0.0015,0.0015,0.0015])
#end_python:

#geometry_objects_read: 0 0 0.00708 C:/Free_CAD/WORKSPACE/Export_Files/vivaldi_antenna_Final/d=0.5mm/antenna_final_d=0.5.h5 C:/AUXILIAR_DOCUMENTOS/Proyecto_Desminado/WORKSPACE/Files/Input_Files/File_Materials/materials_001.txt

#geometry_view: 0 0 0 0.1 0.1 0.05 0.0015 0.0015 0.0015 vivaldi_001_view n