javier/PVPlant
1
0
Fork 0
Code Issues Pull Requests 1 Actions Packages Projects Releases Wiki Activity

updates

Browse Source
This commit is contained in:
Javier Braña
2025-11-20 00:57:15 +01:00
parent 049898c939
commit d61260fdd3
10 changed files with 1704 additions and 520 deletions
Download Patch File Download Diff File Expand all files Collapse all files

214
PVPlantTerrain.py
View File

@@ -73,6 +73,42 @@ line_patterns = {
"Dot (.5x) ...............................": 0x5555,
"Dot (2x) . . . . . . . . . . .": 0x8888}
def open_xyz_mmap(archivo_path):
"""
Usa memory-mapping para archivos muy grandes (máxima velocidad)
"""
# Primera pasada: contar líneas válidas
total_puntos = 0
with open(archivo_path, 'r') as f:
for linea in f:
partes = linea.strip().split()
if len(partes) >= 3:
try:
float(partes[0]);
float(partes[1]);
float(partes[2])
total_puntos += 1
except:
continue
# Segunda pasada: cargar datos
puntos = np.empty((total_puntos, 3))
idx = 0
with open(archivo_path, 'r') as f:
for linea in f:
partes = linea.strip().split()
if len(partes) >= 3:
try:
x, y, z = float(partes[0]), float(partes[1]), float(partes[2])
puntos[idx] = [x, y, z]
idx += 1
except:
continue
return puntos
def makeTerrain(name="Terrain"):
obj = FreeCAD.ActiveDocument.addObject("Part::FeaturePython", "Terrain")
obj.Label = name
@@ -81,7 +117,6 @@ def makeTerrain(name="Terrain"):
FreeCAD.ActiveDocument.recompute()
return obj
class Terrain(ArchComponent.Component):
"A Shadow Terrain Obcject"
@@ -161,101 +196,110 @@ class Terrain(ArchComponent.Component):
if prop == "DEM" or prop == "CuttingBoundary":
from datetime import datetime
if obj.DEM and obj.CuttingBoundary:
'''
Parámetro Descripción Requisitos
NCOLS: Cantidad de columnas de celdas Entero mayor que 0.
NROWS: Cantidad de filas de celdas Entero mayor que 0.
XLLCENTER o XLLCORNER: Coordenada X del origen (por el centro o la esquina inferior izquierda de la celda) Hacer coincidir con el tipo de coordenada y.
YLLCENTER o YLLCORNER: Coordenada Y del origen (por el centro o la esquina inferior izquierda de la celda) Hacer coincidir con el tipo de coordenada x.
CELLSIZE: Tamaño de celda Mayor que 0.
NODATA_VALUE: Los valores de entrada que serán NoData en el ráster de salida Opcional. El valor predeterminado es -9999
'''
grid_space = 1
file = open(obj.DEM, "r")
templist = [line.split() for line in file.readlines()]
file.close()
del file
from pathlib import Path
suffix = Path(obj.DEM).suffix
if suffix == '.asc':
'''
ASC format:
Parámetro Descripción Requisitos
NCOLS: Cantidad de columnas de celdas Entero mayor que 0.
NROWS: Cantidad de filas de celdas Entero mayor que 0.
XLLCENTER o XLLCORNER: Coordenada X del origen (por el centro o la esquina inferior izquierda de la celda) Hacer coincidir con el tipo de coordenada y.
YLLCENTER o YLLCORNER: Coordenada Y del origen (por el centro o la esquina inferior izquierda de la celda) Hacer coincidir con el tipo de coordenada x.
CELLSIZE: Tamaño de celda Mayor que 0.
NODATA_VALUE: Los valores de entrada que serán NoData en el ráster de salida Opcional. El valor predeterminado es -9999
'''
grid_space = 1
file = open(obj.DEM, "r")
templist = [line.split() for line in file.readlines()]
file.close()
del file
# Read meta data:
meta = templist[0:6]
nx = int(meta[0][1]) # NCOLS
ny = int(meta[1][1]) # NROWS
xllref = meta[2][0] # XLLCENTER / XLLCORNER
xllvalue = round(float(meta[2][1]), 3)
yllref = meta[3][0] # YLLCENTER / XLLCORNER
yllvalue = round(float(meta[3][1]), 3)
cellsize = round(float(meta[4][1]), 3) # CELLSIZE
nodata_value = float(meta[5][1]) # NODATA_VALUE
# Read meta data:
meta = templist[0:6]
nx = int(meta[0][1]) # NCOLS
ny = int(meta[1][1]) # NROWS
xllref = meta[2][0] # XLLCENTER / XLLCORNER
xllvalue = round(float(meta[2][1]), 3)
yllref = meta[3][0] # YLLCENTER / XLLCORNER
yllvalue = round(float(meta[3][1]), 3)
cellsize = round(float(meta[4][1]), 3) # CELLSIZE
nodata_value = float(meta[5][1]) # NODATA_VALUE
# set coarse_factor
coarse_factor = max(round(grid_space / cellsize), 1)
# set coarse_factor
coarse_factor = max(round(grid_space / cellsize), 1)
# Get z values
templist = templist[6:(6 + ny)]
templist = [templist[i][0::coarse_factor] for i in np.arange(0, len(templist), coarse_factor)]
datavals = np.array(templist).astype(float)
del templist
# Get z values
templist = templist[6:(6 + ny)]
templist = [templist[i][0::coarse_factor] for i in np.arange(0, len(templist), coarse_factor)]
datavals = np.array(templist).astype(float)
del templist
# create xy coordinates
offset = self.site.Origin
x = (cellsize * np.arange(nx)[0::coarse_factor] + xllvalue) * 1000 - offset.x
y = (cellsize * np.arange(ny)[-1::-1][0::coarse_factor] + yllvalue) * 1000 - offset.y
datavals = datavals * 1000 # Ajuste de altura
# create xy coordinates
offset = self.site.Origin
x = (cellsize * np.arange(nx)[0::coarse_factor] + xllvalue) * 1000 - offset.x
y = (cellsize * np.arange(ny)[-1::-1][0::coarse_factor] + yllvalue) * 1000 - offset.y
datavals = datavals * 1000 # Ajuste de altura
# remove points out of area
# 1. coarse:
if obj.CuttingBoundary:
inc_x = obj.CuttingBoundary.Shape.BoundBox.XLength * 0.0
inc_y = obj.CuttingBoundary.Shape.BoundBox.YLength * 0.0
tmp = np.where(np.logical_and(x >= (obj.CuttingBoundary.Shape.BoundBox.XMin - inc_x),
x <= (obj.CuttingBoundary.Shape.BoundBox.XMax + inc_x)))[0]
x_max = np.ndarray.max(tmp)
x_min = np.ndarray.min(tmp)
# remove points out of area
# 1. coarse:
if obj.CuttingBoundary:
inc_x = obj.CuttingBoundary.Shape.BoundBox.XLength * 0.0
inc_y = obj.CuttingBoundary.Shape.BoundBox.YLength * 0.0
tmp = np.where(np.logical_and(x >= (obj.CuttingBoundary.Shape.BoundBox.XMin - inc_x),
x <= (obj.CuttingBoundary.Shape.BoundBox.XMax + inc_x)))[0]
x_max = np.ndarray.max(tmp)
x_min = np.ndarray.min(tmp)
tmp = np.where(np.logical_and(y >= (obj.CuttingBoundary.Shape.BoundBox.YMin - inc_y),
y <= (obj.CuttingBoundary.Shape.BoundBox.YMax + inc_y)))[0]
y_max = np.ndarray.max(tmp)
y_min = np.ndarray.min(tmp)
del tmp
tmp = np.where(np.logical_and(y >= (obj.CuttingBoundary.Shape.BoundBox.YMin - inc_y),
y <= (obj.CuttingBoundary.Shape.BoundBox.YMax + inc_y)))[0]
y_max = np.ndarray.max(tmp)
y_min = np.ndarray.min(tmp)
del tmp
x = x[x_min:x_max+1]
y = y[y_min:y_max+1]
datavals = datavals[y_min:y_max+1, x_min:x_max+1]
x = x[x_min:x_max+1]
y = y[y_min:y_max+1]
datavals = datavals[y_min:y_max+1, x_min:x_max+1]
# Create mesh - surface:
import MeshTools.Triangulation as Triangulation
import Mesh
stepsize = 75
stepx = math.ceil(nx / stepsize)
stepy = math.ceil(ny / stepsize)
# Create mesh - surface:
import MeshTools.Triangulation as Triangulation
import Mesh
stepsize = 75
stepx = math.ceil(nx / stepsize)
stepy = math.ceil(ny / stepsize)
mesh = Mesh.Mesh()
for indx in range(stepx):
inix = indx * stepsize - 1
finx = min([stepsize * (indx + 1), len(x)-1])
for indy in range(stepy):
iniy = indy * stepsize - 1
finy = min([stepsize * (indy + 1), len(y) - 1])
pts = []
for i in range(inix, finx):
for j in range(iniy, finy):
if datavals[j][i] != nodata_value:
if obj.CuttingBoundary:
if obj.CuttingBoundary.Shape.isInside(FreeCAD.Vector(x[i], y[j], 0), 0, True):
mesh = Mesh.Mesh()
for indx in range(stepx):
inix = indx * stepsize - 1
finx = min([stepsize * (indx + 1), len(x)-1])
for indy in range(stepy):
iniy = indy * stepsize - 1
finy = min([stepsize * (indy + 1), len(y) - 1])
pts = []
for i in range(inix, finx):
for j in range(iniy, finy):
if datavals[j][i] != nodata_value:
if obj.CuttingBoundary:
if obj.CuttingBoundary.Shape.isInside(FreeCAD.Vector(x[i], y[j], 0), 0, True):
pts.append([x[i], y[j], datavals[j][i]])
else:
pts.append([x[i], y[j], datavals[j][i]])
else:
pts.append([x[i], y[j], datavals[j][i]])
if len(pts) > 3:
try:
triangulated = Triangulation.Triangulate(pts)
mesh.addMesh(triangulated)
except TypeError:
print(f"Error al procesar {len(pts)} puntos: {str(e)}")
if len(pts) > 3:
try:
triangulated = Triangulation.Triangulate(pts)
mesh.addMesh(triangulated)
except TypeError:
print(f"Error al procesar {len(pts)} puntos: {str(e)}")
mesh.removeDuplicatedPoints()
mesh.removeFoldsOnSurface()
obj.InitialMesh = mesh.copy()
Mesh.show(mesh)
elif suffix in ['.xyz']:
data = open_xyz_mmap(obj.DEM)
mesh.removeDuplicatedPoints()
mesh.removeFoldsOnSurface()
obj.InitialMesh = mesh.copy()
Mesh.show(mesh)
if prop == "PointsGroup" or prop == "CuttingBoundary":
if obj.PointsGroup and obj.CuttingBoundary: