Files
PVPlant/Utils/PVPlantDownloadOrtophoto.py
2025-03-11 21:18:59 +01:00

97 lines
3.5 KiB
Python

import FreeCAD
import Part
import math
import requests
from io import BytesIO
from PIL import Image
import numpy as np
# ==================================================================
# CONFIGURACIÓN INICIAL (AJUSTAR ESTOS VALORES)
# ==================================================================
lat_min = 40.4165 # Latitud mínima (Ej: Madrid)
lat_max = 40.4300
lon_min = -3.7150 # Longitud mínima
lon_max = -3.6900
zoom = 16 # Nivel de zoom (13-19 para mejor resolución)
terrain_resolution = 50 # Puntos por eje en la malla del terreno
# ==================================================================
def download_tile(lon_deg, lat_deg, zoom):
"""Descarga un tile de OpenStreetMap"""
xtile = int((lon_deg + 180) / 360 * (2 ** zoom))
ytile = int((1 - math.log(math.tan(math.radians(lat_deg)) + 1 / math.cos(math.radians(lat_deg))) / math.pi) / 2 * (
2 ** zoom))
url = f"https://tile.openstreetmap.org/{zoom}/{xtile}/{ytile}.png"
response = requests.get(url)
return Image.open(BytesIO(response.content))
def get_tile_grid(lon_min, lon_max, lat_min, lat_max, zoom):
"""Calcula la cuadrícula de tiles necesarios"""
xtile_min = int((lon_min + 180) / 360 * (2 ** zoom))
xtile_max = int((lon_max + 180) / 360 * (2 ** zoom))
ytile_min = int((1 - math.log(math.tan(math.radians(lat_min)) + 1 / math.cos(math.radians(lat_min))) / math.pi) / 2 * (2 ** zoom))
ytile_max = int((1 - math.log(math.tan(math.radians(lat_max)) + 1 / math.cos(math.radians(lat_max))) / math.pi) / 2 * (2 ** zoom))
return range(xtile_min, xtile_max + 1), range(ytile_min, ytile_max + 1)
def create_stitched_image(xtiles, ytiles, zoom):
"""Une los tiles en una sola imagen"""
tile_size = 256
full_image = Image.new('RGB', (len(xtiles) * tile_size, len(ytiles) * tile_size))
for i, x in enumerate(xtiles):
for j, y in enumerate(ytiles):
tile = download_tile(x, y, zoom)
full_image.paste(tile, (i * tile_size, j * tile_size))
return full_image
def create_parametric_terrain(lon_min, lon_max, lat_min, lat_max, resolution):
"""Crea una superficie paramétrica basada en coordenadas geográficas"""
points = []
for i in np.linspace(lat_min, lat_max, resolution):
row = []
for j in np.linspace(lon_min, lon_max, resolution):
row.append(FreeCAD.Vector(j, i, 0)) # Añadir altura aquí si hay datos de elevación
points.append(row)
terrain = Part.makeFilledFace([points])
obj = FreeCAD.ActiveDocument.addObject("Part::Feature", "Terrain")
obj.Shape = terrain
return obj
def apply_texture(obj, image_path):
"""Aplica la textura al objeto en FreeCAD"""
obj.ViewObject.TextureMode = "EMBEDDED"
obj.ViewObject.TextureImage = image_path
obj.ViewObject.DisplayMode = "Shaded"
# Ejecución principal
'''if __name__ == "__main__":
# Paso 1: Descargar y unir tiles
xtiles, ytiles = get_tile_grid(lon_min, lon_max, lat_min, lat_max, zoom)
stitched_image = create_stitched_image(xtiles, ytiles, zoom)
texture_path = "/tmp/combined_texture.png"
stitched_image.save(texture_path)
# Paso 2: Crear terreno
doc = FreeCAD.newDocument()
terrain = create_parametric_terrain(lon_min, lon_max, lat_min, lat_max, terrain_resolution)
# Paso 3: Aplicar textura
apply_texture(terrain, texture_path)
# Ajustar vista
FreeCAD.Gui.SendMsgToActiveView("ViewFit")
FreeCAD.Gui.updateGui()
print("¡Proceso completado!")'''