helmert-transform 1.2.0

Trasformazione di Helmert 2D per convertire coordinate WGS84 in sistemi locali CAD

Helmert Transform

Libreria Python per la trasformazione di Helmert 2D con shift Z, progettata per convertire coordinate geodetiche WGS84 (longitudine, latitudine, quota) in sistemi di riferimento locali CAD e viceversa.

Caratteristiche

• 🌍 Trasformazione WGS84 ↔ CAD locale bidirezionale usando 5 parametri (TX, TY, TZ, rotazione, scala)
• 🔄 Trasformazione inversa da coordinate locali XYZ a WGS84
• 📐 Ottimizzazione ai minimi quadrati per calcolare i parametri ottimali
• 🎯 Selezione automatica del set di caposaldi più vicino
• 📸 Integrazione Metashape per trasformare coordinate delle foto
• 📊 Statistiche dettagliate sui residui (RMSE, errore massimo, sigma0)
• 🔧 Gestione flessibile di più set di punti di controllo
• 💻 Interfaccia CLI per uso da riga di comando
• 📁 Punti di controllo in CSV per facile gestione e modifica

Installazione

pip install helmert-transform

Da sorgente

git clone https://github.com/marcov/helmert-transform.git
cd helmert-transform
pip install -e .

Utilizzo rapido

Calcolo parametri di trasformazione

from helmert import compute_helmert_parameters, get_control_points

# Ottieni le coordinate dei caposaldi
cp = get_control_points("MU2025")

# Calcola i parametri di trasformazione
result = compute_helmert_parameters(cp['wgs84'], cp['cad'])

# Trasforma un punto WGS84 -> CAD
x, y, z = result['transform'](11.27, 45.59, 500.0)
print(f"Coordinate CAD: ({x:.2f}, {y:.2f}, {z:.2f})")

# Trasformazione inversa: CAD -> WGS84
lon, lat, h = result['inverse_transform'](x, y, z)
print(f"Coordinate WGS84: ({lon:.8f}°, {lat:.8f}°, {h:.2f} m)")

Trasformazione inversa: da coordinate locali a WGS84

from helmert import compute_helmert_parameters, get_control_points

cp = get_control_points("MU2025")
result = compute_helmert_parameters(cp['wgs84'], cp['cad'], verbose=False)

# Lista di punti in coordinate locali CAD
local_points = [
    (1677100.0, 5051000.0, 450.0),
    (1677200.0, 5051100.0, 455.0),
    (1677150.0, 5051050.0, 460.0),
]

# Converti in WGS84
for x, y, z in local_points:
    lon, lat, h = result['inverse_transform'](x, y, z)
    print(f"({x:.2f}, {y:.2f}, {z:.2f}) -> ({lon:.8f}°, {lat:.8f}°, {h:.2f} m)")

Definire punti di controllo personalizzati

from helmert import compute_helmert_parameters

# Coordinate WGS84 (lon, lat, h) dei punti noti
wgs84_points = [
    (11.266, 45.592, 434.0),
    (11.267, 45.591, 442.0),
    (11.269, 45.593, 450.0),
]

# Coordinate corrispondenti nel sistema CAD locale
cad_points = [
    (1676805.01, 5051329.05, 386.22),
    (1676819.70, 5051230.11, 394.72),
    (1676912.45, 5051462.58, 402.09),
]

# Calcola la trasformazione
result = compute_helmert_parameters(wgs84_points, cad_points)

# Statistiche
print(f"RMSE: {result['statistics']['rmse']*1000:.3f} mm")
print(f"Scala: {result['parameters']['scale']:.10f}")
print(f"Rotazione: {result['parameters']['rotation_deg']:.6f}°")

Selezione automatica del set di caposaldi

from helmert import find_best_control_point_set, CONTROL_POINT_SETS

# Coordinate delle immagini (es. da GPS delle foto)
image_coords = [
    (11.270, 45.593, 500.0),
    (11.272, 45.594, 510.0),
    (11.271, 45.592, 495.0),
]

# Trova il set più vicino
result = find_best_control_point_set(image_coords, CONTROL_POINT_SETS)
print(f"Set consigliato: {result['best_set']}")
print(f"Distanza: {result['distance']/1000:.2f} km")

Utilizzo con Agisoft Metashape

# Dalla console Python di Metashape
import sys
sys.path.insert(0, '/path/to/helmert-transform')

from helmert.metashape import main

# Esegui con selezione automatica del set
main(auto_select=True)

# Oppure specifica il set manualmente
main(set_name="MU2025")

# Usa un file CSV esterno con punti di controllo personalizzati
main(control_points_file="/path/to/my_control_points.csv", set_name="MySet")

# Combina file esterno con selezione automatica
main(control_points_file="/path/to/my_control_points.csv", auto_select=True)

Interfaccia da riga di comando (CLI)

Il package include un'interfaccia da riga di comando per operazioni comuni.

Installazione CLI

Dopo l'installazione del package, il comando helmert sarà disponibile:

pip install -e .
helmert --help

Opzione globale: file CSV esterno

È possibile utilizzare un file CSV esterno con punti di controllo personalizzati, invece del file predefinito incluso nel package:

# Usa un file CSV esterno
helmert --control-points /path/to/my_control_points.csv list-sets

# Abbreviazione: --cp
helmert --cp my_points.csv transform --set MySet --lon 11.27 --lat 45.59 --h 500

Il file CSV esterno deve avere il formato:

set_name,point_id,lon,lat,h,x,y,z
MySet,P01,11.266,45.592,434.0,1676805.01,5051329.05,386.22
MySet,P02,11.267,45.591,442.0,1676819.70,5051230.11,394.72

Comandi disponibili

Elencare i set di caposaldi

helmert list-sets

# Con file CSV esterno
helmert --cp my_points.csv list-sets

Informazioni su un set

helmert info --set MU2025

# Con file CSV esterno
helmert --cp my_points.csv info --set MySet

Mostra i dettagli del set, le coordinate e i parametri di trasformazione calcolati.

Trasformare un singolo punto (WGS84 → CAD)

# Output compatto
helmert transform --lon 11.27 --lat 45.59 --h 500

# Output dettagliato
helmert transform --set MU2025 --lon 11.27 --lat 45.59 --h 500 --verbose

Trasformazione inversa (CAD → WGS84)

# Singolo punto
helmert inverse --x 1677100 --y 5051000 --z 450

# Con output dettagliato
helmert inverse --set MU2025 --x 1677100 --y 5051000 --z 450 --verbose

Trasformazione batch da file (WGS84 → CAD)

# Da file a file
helmert batch --set MU2025 --input punti_wgs84.txt --output punti_cad.txt

# Da stdin a stdout (per pipeline)
cat punti.txt | helmert batch --set MU2025

Formato del file di input (un punto per riga, ID come prima colonna):

id lon lat h
P001 11.270 45.593 500.0
P002 11.272 45.594 510.0

Formato output:

id x y z
P001 1677102.566 5051032.032 451.953
P002 1677250.123 5051145.678 461.953

Trasformazione batch inversa (CAD → WGS84)

# Da file a file
helmert batch-inverse --set MU2025 --input punti_cad.txt --output punti_wgs84.txt

Formato del file di input (un punto per riga, ID come prima colonna):

id x y z
P001 1677100.0 5051000.0 450.0
P002 1677200.0 5051100.0 455.0

Formato output:

id lon lat h
P001 11.26994275 45.58971368 498.047
P002 11.27129833 45.59055829 503.047

Esportazione in KML

# Converte punti CAD in file KML per Google Earth
helmert to-kml --input punti_cad.txt --output punti.kml

# Con nome personalizzato
helmert to-kml -i punti.txt -o output.kml --name "Rilievo 2025"

Il formato di input è lo stesso di batch-inverse (id x y z).

Esportazione control points di un set

# Esporta le coordinate WGS84 dei caposaldi di un set in KML
helmert export-set --set MU2025 --output caposaldi.kml

# Con nome personalizzato per il documento
helmert export-set -s MU2025 -o caposaldi.kml --name "Caposaldi Rilievo 2025"

Questo comando è utile per visualizzare i punti di controllo in Google Earth.

Struttura del package

helmert-transform/
├── helmert/
│   ├── __init__.py         # Esportazioni pubbliche
│   ├── core.py             # Funzioni di calcolo Helmert
│   ├── control_points.py   # Gestione punti di controllo
│   ├── cli.py              # Interfaccia riga di comando
│   └── metashape.py        # Integrazione Agisoft Metashape
├── data/
│   └── control_points.csv  # Punti di controllo in formato CSV
├── tests/                  # Test unitari
├── pyproject.toml          # Configurazione package
└── README.md

API Reference

Funzioni principali

compute_helmert_parameters(wgs84_points, cad_points, verbose=True)

Calcola i parametri ottimali della trasformazione di Helmert.

Parametri:

• wgs84_points: Lista di tuple (lon, lat, h) in gradi e metri
• cad_points: Lista di tuple (X, Y, Z) in metri
• verbose: Se True, stampa informazioni dettagliate

Ritorna: Dizionario con:

• parameters: TX, TY, TZ, rotazione, scala
• statistics: RMSE, errori, residui
• transform: Funzione di trasformazione diretta (lon, lat, h) -> (X, Y, Z)
• inverse_transform: Funzione di trasformazione inversa (X, Y, Z) -> (lon, lat, h)

get_control_points(set_name=None)

Restituisce un set di caposaldi.

Parametri:

• set_name: Nome del set (default: set attivo)

Ritorna: Dizionario con point_ids, wgs84, cad

find_best_control_point_set(image_coords, control_point_sets, verbose=True)

Trova il set di caposaldi più vicino geograficamente.

Parametri:

• image_coords: Lista di coordinate WGS84 delle immagini
• control_point_sets: Dizionario dei set disponibili
• verbose: Se True, stampa informazioni

Ritorna: Dizionario con best_set, distance, bbox

Costanti

• WGS84_A: Semi-asse maggiore WGS84 (6378137.0 m)
• WGS84_F: Appiattimento WGS84 (1/298.257223563)
• WGS84_E2: Eccentricità al quadrato

Come aggiungere nuovi set di caposaldi

I punti di controllo sono definiti nel file data/control_points.csv. Il formato è:

set_name,point_id,lon,lat,h,x,y,z
MU2025,M01,11.26641983,45.59282693,434.268,1676805.01,5051329.05,386.22
MU2025,M02,11.26653239,45.59193042,442.774,1676819.70,5051230.11,394.72
...
nuovo_set,P01,12.345,46.789,500.0,1234567.89,9876543.21,450.0
nuovo_set,P02,12.346,46.790,510.0,1234600.00,9876600.00,460.0

Colonne del file CSV

Colonna	Descrizione
set_name	Nome del set di caposaldi (es. MU2025, ridotto)
point_id	Identificativo del punto (es. M01, P001)
lon	Longitudine WGS84 in gradi decimali
lat	Latitudine WGS84 in gradi decimali
h	Quota ellissoidica WGS84 in metri
x	Coordinata X nel sistema CAD locale (metri)
y	Coordinata Y nel sistema CAD locale (metri)
z	Coordinata Z nel sistema CAD locale (metri)

Aggiungere un nuovo set

Per aggiungere un nuovo set di caposaldi:

1. Apri il file data/control_points.csv
2. Aggiungi le righe con il nuovo nome del set nella prima colonna
3. Servono almeno 3 punti per calcolare la trasformazione

Le righe che iniziano con # sono considerate commenti.

Caricare set da file esterno

È possibile caricare punti di controllo da un file CSV esterno:

from helmert.control_points import load_and_register_csv

# Carica e registra i set dal file
load_and_register_csv('/path/to/mio_file.csv')

# Ora i set sono disponibili
from helmert import get_control_points
cp = get_control_points("nome_set_nel_csv")

Teoria

La trasformazione di Helmert 2D (similitudine piana) con shift Z utilizza 5 parametri:

1. TX, TY: Traslazioni nel piano
2. TZ: Traslazione verticale
3. θ: Rotazione nel piano
4. k: Fattore di scala

La trasformazione avviene in due passi:

1. WGS84 → ENU: Le coordinate geodetiche vengono proiettate su un piano tangente locale (East-North-Up) centrato sul centroide dei punti
2. ENU → CAD: Applicazione della trasformazione di Helmert

I parametri vengono calcolati minimizzando la somma dei quadrati dei residui usando l'algoritmo di Levenberg-Marquardt.

Perché Helmert 2D a 5 parametri e non 3D a 7?

Una trasformazione di Helmert 3D completa utilizza 7 parametri:

• 3 traslazioni (TX, TY, TZ)
• 3 rotazioni (ωx, ωy, ωz) attorno ai tre assi
• 1 fattore di scala

Questa libreria implementa invece una Helmert 2D con shift Z a 5 parametri per le seguenti ragioni:

1. Assi Z quasi paralleli: Nel caso d'uso tipico (rilievo fotogrammetrico da drone), l'asse Z del sistema WGS84 (verticale ellissoidica) e l'asse Z del sistema CAD locale (verticale locale) sono praticamente paralleli. La differenza è data solo dalla deflessione della verticale, che in Italia è nell'ordine di pochi secondi d'arco (trascurabile per applicazioni topografiche).

2. Minori gradi di libertà: Con 5 parametri invece di 7, servono meno punti di controllo per una soluzione stabile. Bastano 3 punti ben distribuiti, mentre per 7 parametri ne servirebbero almeno 4.

3. Stabilità numerica: Meno parametri significano meno correlazioni tra le incognite e una stima più robusta, specialmente quando i punti di controllo coprono un'area limitata.

4. Coerenza con la pratica topografica: La maggior parte dei software CAD e GIS per rilievi locali utilizza trasformazioni piane (similitudine o affine) con gestione separata della quota.

5. Adattamento altimetrico: Il parametro TZ assorbe la differenza sistematica tra quota ellissoidica WGS84 e quota del sistema locale (che include l'ondulazione del geoide, tipicamente ~48m nel Nord Italia).

Per trasformazioni tra datum geocentrici diversi (es. WGS84 ↔ ED50) sarebbe invece necessaria una Helmert 3D completa a 7 parametri.

Requisiti

• Python ≥ 3.8
• NumPy ≥ 1.20.0
• SciPy ≥ 1.7.0
• (Opzionale) Agisoft Metashape per l'integrazione fotogrammetrica

Licenza

Questo progetto è distribuito sotto licenza MIT. Vedi LICENSE per i dettagli.

Contribuire

1. Fork del repository
2. Crea un branch per la feature (git checkout -b feature/nuova-feature)
3. Commit delle modifiche (git commit -am 'Aggiunge nuova feature')
4. Push del branch (git push origin feature/nuova-feature)
5. Apri una Pull Request

Changelog

Vedi CHANGELOG.md per la storia delle versioni.