Tutorial esercitazioni inquinamento elettromagnetico a bassa frequenza Aggiornato all' 11 marzo 2008

Indice

• Introduzione
• Link utili
• Impostazioni preliminari
  • Generazione proiezione GAUSS-BOAGA
  • Proprietà del progetto
  • Generazione nuovo Database / Location / Mapset
• Visualizzazione Raster, vista panoramica e bookmarks
• Modello del terreno
  • Impostazione della regione di lavoro di GRASS
  • Importazione shapefile punti quotati e curve di livello.
  • Trasformazione degli shapefile in in layer vettoriali GRASS.
  • Trasformazione di un layer vettoriale di linee in punti e unione di due layer.
  • Generazione di un DTM Raster a partire da un layer vettoriale GRASS di punti.
  • Modifica della mappa di colori del raster GRASS che rappresenta il DTM.
  • Terrain analysis.
• Digitalizzazione dell'elettrodotto
  • Generazione di un nuovo shapefile puntuale dei sostegni
  • Importazione e visualizzazione sostegni in GRASS (etichette).
  • Determinazione quota dei sostegni in base a DTM e aggiunta coordinate X,Y,Z sostegni nella tabella degli attributi.
  • Generazione tracciato (vettore lineare) da sostegni (vettore di punti): approccio manuale
  • Generazione tracciato (vettore lineare) da sostegni (vettore di punti): approccio automatizzato (con script v.to.lines)
• Importazione e primo filtraggio edificato
• Generazione di nuovi location / mapset e riproiezione di layer in diversi sistemi di coordinate
• Generazione nuovi Location / Mapset nello stesso Database
• Riproiezione di un vettore
• Determinazione buffer e grigliato
  • Buffer intorno al tracciato della linea
  • Grigliato regolare di punti
  • Estrazione dei punti del grigliato all'interno del buffer a 300 m
  • Quotatura grigliato in base a DTM
• Modellazione elettrodotto e calcolo del campo magnetico
  • Modellazione elettrodotto
  • Calcolo del campo magnetico
• Modellazione elettrodotto e calcolo del campo magnetico
  • Visualizzazione shp puntuale con "simbolo graduato", "colore continuo", "valore univoco"
  • Individuazione edifici con altezza superiore a 7 m con livelli di esposizione localmente maggiori di 3 microtesla
  • Suddivisione dell'edificato con altezza superiore a 7 m in porzioni di edificio differenziate in base ai livelli di esposizione, approccio 1
  • Suddivisione dell'edificato con altezza superiore a 7 m in porzioni di edificio differenziate in base ai livelli di esposizione, approccio 2

Questo tutorial è stato creato come supporto alle 10 ore di esercitazione tenute da Nicola Zoppetti nell'ambito del Master di secondo livello "Sistemi Informativi Geografici per il monitoraggio e la gestione del territorio" 2006-2007. Le esercitazioni seguono 4 ore di seminario sul tema dell'inquinamento elettromagnetico a bassa frequenza ed in particolare sulle linee elettriche ad alta tensione. Il tutorial è stato ulteriormente sviluppato nell'ambito dell'edizione 2008 del master. Nel seminario ci si è concentrati sugli aspetti fisici e tecnici riguardanti gli elettrodotti ed i campi da essi generati. Nelle esercitazioni invece si affronta un caso studio utilizzando strumenti GIS. L'aspetto del calcolo dei campi viene solo accennato, trattandosi di un tema che richiederebbe ben più di 10 ore per essere affrontato. Ci si concentra invece sulla preparazione dei dati necessari al calcolo e sull'elaborazione dei risultati. Per vari motivi si è scelto di sviluppare le esercitazioni vincolandosi al sistema operativo Microsoft Windows.
Il carattere della esercitazione vuole essere estremamente applicativo: non ci si soffermerà su aspetti teorici o tecnici, pur fondamentali, riguardanti i GIS in generale o QGIS-GRASS in partcolare. Piuttosto si illustrerà come possono essere raggiunti determinato obiettivi rimandando per gli approfondimenti e per gli aspetti generali o di principio, ad opportuni riferimenti.

Per le esercitazioni si è scelto di utilizzare il software QGIS. Si tratta di un software Opne Source, disponibile gratuitamente per vari sistemi operativi. Il software può essere scaricato qui. Il programma è stato originariamente sviluppato per linux e solo recentemente è possibile utilizzarlo anche in windows. La versione per windows, essendo meno matura è anche meno stabile di quella per linux e quindi ci sono alcune funzionalità ancora non supportate o che non funzionano come dovrebbero. Si tratta in ogni caso di un programma in continua evoluzione e dalle grandi potenzialità e quindi vale la pena provare ad utilizzarlo anche se non si può dire perfettamente funzionante. Attualmente QGIS, che è nato come GIS leggero e con funzionalità di base, può essere utilizzato insieme a GRASS, un GIS tra i più complessi ed evoluti, e ciò ha ampliato di molto le possibilità di analisi. Sotto windows, quanto detto per QGIS, vale a maggior ragione anche per GRASS, visto che il suo "porting" per windows è ancora più recente. Al momento, ad esempio, non è supportata la visualzzazione tridimensionale.
Questo tutorial verrà auspicabilmente aggiornato ed esteso sia per illustrare eventuali nuove funzionalità, sia per affrontare i molti aspetti e le possibili elaborazioni legate al problema della valutazione dell'impatto ambientale del campo elettrico e magnetico generati da elettrodotti che non si è avuto tempo di trattare in questa prima versione.

• Sito ufficiale QGIS
• Blog QGIS
• Download ultima versione QGIS
• Wiki QGIS: indice dei titoli
• Archivi mailing-list sviluppatori di QGIS (per controllare lo stato di sviluppo)
• Archivi mailing-list utenti di QGIS

• GFOSS ITALIA: Geospatial Free Open Source Software
• GFOSS: HOWTOs e procedure tipo

• Sito ufficiale GRASS
• Sito utenti GRASS italia
• Lista utenti GRASS italia
• Ottimo tutorial GRASS
• Manuali varie versioni di GRASS

• GRASS per windows
• Download WINGRASS
• WinGRASS Current Status

Installare QGIS ( QGIS ) in una cartella il cui percorso sia privo di spazi (ad es. C:/Programmi/QGIS)

Creare una nuova cartella di lavoro per le esercitazioni. Ad esempio c:/esercitazioni/

Se si vuole salvare il tutorial , scaricare il file Tutorial HTML (zip 9 MB) ed estrarre tutti i files nella cartella appena creata (utilizzare l'opzione "Use folder names" di WinZip) Nota: i link ai file .zip non funzioneranno nella versione locale.

Scaricare il file Dati.zip (zip 24 MB) ed estrarre tutti i files nella cartella esercitazioni (utilizzare l'opzione "Use folder names" di WinZip)

GRASS, come molti altri sistemi GFOSS, usa la libreria PROJ.4 per gestire i sistemi di riferimento, le proiezioni e le loro trasformazioni. Le PROJ.4 definiscono i sistemi di riferimento in base alle loro differenze rispetto al WGS84 (“il datum del GPS”): impostare un sistema di riferimento significa quindi scegliere i parametri della trasformazione da un sistema al WGS84. In teoria si usa una trasformazione a 7 parametri (tre traslazioni + tre rotazioni + 1 fattore di scala), per come sono realizzati i sistemi di riferimento spesso si hanno errori inaccettabili e si usano griglie di valori, i parametri cambiano cioè da punto a punto. link www alla presentazione pdf "Sistemi di riferimento e proiezioni in GRASS" di Paolo Zatelli

Dopo aver lanciato per la prima volta QGIS, impostare un nuovo sistema di proiezione: menù Impostazioni / Proiezione personalizzata / scheda Definizione

Nome: Gauss-Boaga corretto Parametri: +proj=tmerc +lat_0=0 +lon_0=9 +k=0.999600 +x_0=1500000 +y_0=0 +ellps=intl +units=m +towgs84=-104.1,-49.1,-9.9,0.971,-2.917,0.714,-11.68

menù Impostazioni / Proprietà del progetto /

scheda Generale unità mappa: metri precisione: manuale, 3 posizioni decimali

menù Impostazioni / Proprietà del progetto / scheda Proiezione abilita la proiezione al volo (check box in alto) user defined ... / Gauss-Boaga corretto

Salvare il progetto menù file / salva

menù Plugin / Grass / New Mapset

Creare una nuova cartella "Datigrass" (ad esempio in c:/esercitazioni/Dati) Next

Creare una nuova location "ToscanaGaussBoaga" Next

Selezionare la proiezione definita precedentemente Next

Impostare la regione GRASS sull'Italia Next

Nuovo mapset :"nicola" Next

Finish

menù Layer / aggiungi un raster

selezionare 275020.tif 275030.tif

menù Visualizza / vista massima

menù Layer / aggiungi tutto alla vista panoramica

individuare la zona di interesse (vedere riquadro nella vista panoramica)

menù Visualizza / nuovo segnalibro : "Argingrosso"

D'ora in poi, ogni volta che si vuole visualizzare la zona di interesse sarà sufficiente utilizzare il segnalibro appena creato: menù Visualizza / Mostra segnalibri...

Zoom a

La regione di interesse è coperta dal solo layer 275030 per cui si può eliminare il layer 275020. Eliminare non vuol dire cancellare il file dal disco ma semplicemente non prenderlo in considerazione nel progetto corrente. selezionare il layer nella legenda: menù Layer / elimina layer oppure ctrl+D

menù Plugins / Grass / apri strumenti Grass

Toolbox GRASS /Modules /File / Import / Import raster /r.in.GDAL

Esegui

Chiudi

Toolbox GRASS / Browser / Imposta la regione corrente con la mappa selezionata

Selezionare 275030 nella legenda e Zoom al layer

Selezionare 275030 (il .tif) nella legenda e Elimina Layer

Aggiungi raster Grass / 275030

Selezionare 275030 nella legenda e Nella vista panoramica

Salvare il progetto menù file / salva

caricare i layer shp vettoriali delle curve di livello e dei punti quotati menù Layer / Aggiungi layer vettoriale

CARTOGRAFIA_BASE/275030_PuntiQuotati.shp CARTOGRAFIA_BASE/275030_CurveLivello.shp

Cambiare ordine di visualizzazione nella legenda: trascinare con il puntatore del mouse il raster 275030 in cima alla lista (scompaiono i punti e le curve di livello)

Doppio clic nella legenda su 275030 (oppure selezionare 275030 e poi Tasto destro / Proprietà) Trasparenza 50% OK

selezionare nella legenda 275030_CurveLivello Toolbox GRASS /Modules /File / Import / Import vector /v.in.ogr

OGR vector layer: 275030_CurveLivello.shp name for output...:c275030

Esegui

OGR vector layer: 275030_PuntiQuotati.shp name for output...:p275030 Esegui

Selezionare gli shp nella legenda (uno alla volta) layer /elimina layer

aggiungere i layer vettoriali GRASS al prgetto ci sono 3 modi: Plugins / GRASS / aggiungi vettoriale GRASS icona aggiungi vettoriale GRASS Toolbox GRASS / Browser / icona Aggiungi la mappa selezionata all'area di lavoro

Toolbox GRASS /Modules /File / Map type conversion / convert a raster to vector / v.to.points Mella versione 0.9.2 Toolbox GRASS /Modules / vector / geometry management /manage features / v.to.points

input: c275030 output:pc275030 distance 10 m Esegui

Toolbox GRASS / Browser / icona Aggiungi la mappa selezionata all'area di lavoro pc275030 1 point

eliminare layer c275030

unione di due layer puntuali: primo approccio Toolbox GRASS / Modules /GRASS shell Compare una shell (finestra comandi) scrivere nella shell: v.patch input=pc275030,p275030 output=ppc275030_1 copiare la riga sopra e incollarla nella shell (maiusc+tasto sinistro mouse per incollare nella finestra dei comandi GRASS)

Il risultato è una mappa con due layer (perché?) composta dall' unione dei punti delle due mappe di partenza. La tabella degli attributi non ha conservato i campi delle tabelle originarie e quindi ha perso anche il riferimento alla quota.

Si prova ad aggiungere al comando l'opzione -e che indica di copiare anche la tabella degli attributi. Inoltre specifico anche l'opzione --overwrite in modo da poter lanciare il comando ripetutamente utilizzando la memoria della finestra dei comandi. v.patch input=pc275030,p275030 output=ppc275030_2 -e --overwrite

Si prova ad aggiungere al comando l'opzione -e che indica di copiare anche la tabella degli attributi. Inoltre specifico anche l'opzione --overwrite in modo da poter lanciare il comando ripetutamente utilizzando la memoria della finestra dei comandi. v.patch input=pc275030,p275030 output=ppc275030_2 -e --overwrite

Il risultato dell'esecuzione di questo comando evidenzia come GRASS sotto windows non sia ancora perfettamente funzionante. A volte il comando fallisce:

A volte il comando va a buon fine:

Quando il comando va a buon fine si ottiene una mappa con un solo layer puntuale a cui è associata una tabella degli attributi completa.

Unione di due layer puntuali: approccio alternativo trasformare i layer dei punti GRASS in SHP (problemi, si pianta se uso v.out.ogr) tasto destro nella lista/salva come shp importare in un unico layer tutti i punti: creare una nuova cartella (ad esempio c:/archivio/mastergis/esercitazioni/work/punti) v.in.ogr dsn=c:/archivio/mastergis/esercitazioni/work/punti output=puntiquotati

visualizzare la regione corrente (riquadro rosso)

icona modifica regione grass attuale o del menù plugins/grass/modifica regione grass attuale) impostare la regione di lavoro adattata al layer dei punti considerato impostare una risoluzione di 10m sia per le righe sia per le colonne. La risoluzione ed il corrispondente numero di righe e colonne nella regione di lavoro, determinano la dimensione dei raster generati con gli strumenti di grass. Se si scelgono risoluzioni che comportano immagini con molti pixel, i tempi di calcolo si allungano di molto ed aumentano le risorse di memoria necessarie). Con normali pc si hanno tempi di calcolo ragionevoli con raster di 1000x1000 pixel.

Toolbox GRASS / Modules / raster / interpolate surface v.surf.idw nella versione 0.9.0 (a differenza che nella 0.8.1) non è possibile scegliere mediante menù a discesa l'attributo in base al quale generare il raster dtm. Nella nuova versione 0.9.2rc1 ciò è di nuovo possibile.

Si ricorre ancora una volta alla shell di GRASS Toolbox GRASS / Modules /GRASS shell v.surf.idw input=ppc275030_3 column="VALORE" output=dtm_10m_8p npoints=8

Aggiungere il raster grass creato (icona o menù Plugin / GRASS /aggiungi raster)

Si crea un nuovo dtm utilizzando un numero diverso di punti (4 in questo caso) nella GRASS shell v.surf.idw input=ppc275030_3 column="VALORE" output=dtm_10m_4p npoints=4

Si crea un nuovo raster in cui sono riportate l differenze tra i due dtm: nella GRASS shell r.mapcalc 'delta_dtm_idw8pt_idw4pt=dtm_10m_8p-dtm_10m_4p'

Visualizzare l'istogramma di uno dei dtm, ad esempio dtm_10m_8p. Tasto destro / proprietà / istogramma / tasto aggiorna. Si nota che le quote sono comprese tra 20 e 160 m s.l.m. provare a cambiare ilò numero delle colonne

Utilizzare il comando r.colors (lanciarlo da shell) utilizzando "valori percentuali" r.colors map=dtm_10m_8p color=rules 0% 0 0 0 20% 0 100 0 40% 0 255 0 60% 0 100 100 80% 0 255 255 100% 255 0 0 end

Invece di scrivere sulla shell le regole, è possibile memorizzarle una volta per tutte su un file di testo In particolare creare un file di nome color1.file con il seguente contenuto 0% 0 0 0 20% 0 100 0 40% 0 255 0 60% 0 100 100 80% 0 255 255 100% 255 0 0 end

Da shell, lanciare il comando: cat C:/Archivio/MasterGIS/Esercitazioni/work/color1.file | r.colors map=dtm_10m_4p color=rules

E' possibile lanciare il comando dall'interfaccia grafica del gestore strumenti di GRASS: in tal caso eliminare end finale dal file delle regole

E' inoltre possibile utilizzare il comando r.colors utilizzando "valori assoluti" r.colors map=dtm_10m_8p color=rules 0 0 0 255 50 0 255 0 70 0 100 0 100 50 100 0 130 100 100 0 200 255 100 0 end

Ci si può poi divertire utilizzando i vari strumenti GRASS dedicati alla "Terrain Analysis". ad esempio è possibile generare la cosiddetta mappa shaded cioè una visualizzazione del modello del terreno in cui compaiono le ombre che ne evidenziano l'orografia.

Impostando opportunamente la trasparenza dei raster (doppio clic sul nome nella legenda o tasto destro + proprietà) è possibile sovrapporre più immagini e ottenere mappe composite (nell'immagine la sovrapposizione tra il dem, la mappa shaded, quella delle pendenze e quella degli "aspect").

Scaricare o spengere tutte le mappe. Caricare il raster 275030.tif (non il raster GRASS): menù layer / aggiungi un raster Doppio clic sul nome della mappa nella legenda oppure tasto destro / proprietà scheda piramidi evidenziare ad una ad una le risoluzioni desiderate e premere crea piramidi (a volte il programma va in crash ma a forza di tentativi ci si riesce)

menù Layer / nuovo layer vettoriale

tipo: punto formato: ESRI shapefile attributi: X real Y real Z real h real tipo string codice string

OK salvare lo shapefile nella cartella work con nome sostegni_407

visualizza / mostra segnalibri / Argingrosso / zoom a seguire il tracciato dell'elettrodotto prima verso est e poi verso sud fino alla sottostazione presso "il ronco".

icona toggle editing icona inserisci punti con il mouse cliccare nel punto in cui sul raster è rappresentato il primo sostegno della linea (all'interno della sottostazione di S.Lorenzo a Greve) inserire gli attributi nella tabella (lasciare vuoti i campi x,y,z).

modificare la simbologia dell' shp in modo che sia più evidente.

utilizzando una scala di visualizzazione che permetta l'inserimento dei sostegni sulla mappa con una precisione sufficiente, alternare l'inserimento dei nuovi sostegni e il pan lungo la linea.

una volta arrivati al fiume arno terminare l'inserimento icona toggle editing salvare le modifiche

evidenziare sostegni_407 nella legenda menù Visualizza / zoom sul layer

evidenziare sostegni_407 nella legenda tasto destro / visualizza tabella attributi oppure icona Apri tabella

correggere gli eventuali errori di inserimento Abilita modifica doppio clic sul campo da modificare alla fine di tutte le modifiche Interrompi modifica Salva Apply Ok

importare l'shp sostegni_407 come layer GRASS apri strumenti GRASS / v.in.ogr eliminare il layer sostegni_407.shp caricare il layer grass sostegni_407 (icona o menù Plugin / GRASS / Aggiungi vettoriale GRASS)

Aggiunta etichette ai sostegni evidenziare sostegni_407 nella legenda tasto destro / proprietà oppure doppio clic sul nome del layer Scheda etichette: sub-scheda posizione: campo con etichetta: codice_tar spostamento unità: punti sub-scheda carattere: impostare a piacere sub-scheda allineamento carattere: a destra nell'esempio in figura

caricare il raster GRASS dtm_090_10m creato in precedenza

Toolbox GRASS Modules / Vectors / Develop map / convert 2D vector to 3D by sampling raster

input: sostegni_407 output: sostegni_407_3d raster: dtm_090_10m sampling: a piacere, se c'è tempo fare un po' di prove

Caricare il layer vettoriale GRASS sostegni_407_3d Toolbox GRASS Modules / Vectors / Report and statistics / v.to.db

Con questo comando è possibile agiungere alla tabella degli attributi le coordinate dei punti. l'interfaccia del comando v.to.db non permette di specificare i campi in cui inserire le coordinate (non funziona il menù a discesa).

si usa la shell v.to.db map=sostegni_407_3d option=coor units=me column=X,Y,Z

NB: il baco descritto di seguito è stato risolto nell'ultima versione probabilmente a causa di un "baco" il programma non aggiorna la tabella dbf associata al layer vettoriale ma crea solo una tabella temporanea. si deve sostituire manualmente la tabella associata con quella temporanea. Per farlo: si naviga con esplora risorse fino alla cartella .tmp nel mapset in cui stiamo lavorando. Si ordinano i file secondo la data di ultima modifica Si copia in clipboard l'ultimo dbf creato Si incolla il dbf nella cartella dbf, allo stesso livello di .tmp si rinomina: sostegni_407_3d ---> sostegni_407_3d_old si rinomina 2456 (file temporaneo) ---> sostegni_407_3d Si visualizza la tabella degli attributi in cui i campi X, Y e Z sono correttamente riempiti (Se ciò non accade provare a scaricare e ricaricare il layer sostegni_407_3d affinché l'applicazione si "accorga" delle modifiche alla tabella dbf.

tasto destro nella legenda su sostegni_407_3d esporta come shapefile salvare sostegni_407_3d nella cartella work un'ulteriore opportunità di esportazione è data dallo strumento grass toolbox GRASS / File / Export / Export vectors / v.out.ascii

lanciando il comando da shell è possibile personalizzare il carattere separatore ed altri particolari. salvare il layer grass sostegni_407_3d in file di testo sia in formato point sia in formato standard salvare in due file diversi nella cartella work

aprire con wordpad il file s407_point.txt sostitire il carattere | con uno spazio per evitare questo passaggio si poteva lanciare da shell il comando di esportazione in formato testo specificando il separatore spazio o tab: v.out.ascii input=sostegni_407_3d format=point output=C:/Archivio/MasterGIS/Esercitazioni/Dati/DatiElettrodotti/sostegni_407_3d.txt fs="\t"

aprire con wordpad il file s407_std.txt cancellare tutto il contenuto del file dopo la 10 riga VERTI(esclusa) aggiungere la riga L 16 1 copiare dopo l'ultima riga tutto il contenuto di s407_point.txt inserire come ultima riga: 1 1

dalla shell di GRASS lanciare il comando v.in.ascii input=C:/Archivio/MasterGIS/Esercitazioni/Dati/work/s407_std.txt output=tracciato407 --overwrite format=standard

caricare il layer vettoriale tracciato407

Il precedente approccio può essere automatizzato mediante la generazione di uno script. Non si sta qui ad entrare nel merito di come è possibile generare uno script da utilizzare in grass. Lo script già fatto si trova nella cartella Dati/varie e si chiama v.to.lines Copiare il file ed incollarlo nella cartella ..../QuantumGIS/grass/scripts Dove ..../QuantumGIS è la cartella di installazine del programma Una volta copiato lo script provare a digitare nella shell: v.to.lines

Lo script richiede che si specifichi: il nome del layer puntuale di ingresso il nome del layer delle linee in uscita

v.to.lines input=sostegni_407_3d output=t407

Per generare lo script v.to.lines ho utilizzato i seguenti link: Messagio 1 nella lista utenti grass Messagio 2 nella lista utenti grass Tutorial awk Quick card reference awk

Se poi si vuole aggiungere ala toolbox grass l'interfaccia grafica riferirsi al seguente link: QGIS: Grass Tools aggiunta moduli

In particolare: copiare i file v.to.lines.qgs v.to.lines.1.png e v.to.lines.2.png nella cartella ..../QuantumGIS/grass/modules modificare il file ..../QuantumGIS/grass/config/default.qgc agiungendo sotto grass name="v.to.points"/ la riga grass name="v.to.lines"/

Affinché le modifiche abbiano effetto, uscire da QGIS e rilanciarlo

menù Plugin / Grass / New Mapset

Riferirsi allo stesso Database "Datigrass" (e quindi cartella) creato in precedenza Next

Creare una nuova location "ToscanaUTMED50" Next

Selezionare la proiezione con id 1250 Next

Nuovo mapset :"nicola" Next

menù Plugin / Grass / New Mapset

Riferirsi allo stesso Database "Datigrass" (e quindi cartella) creato in precedenza Next

Creare una nuova location "ToscanaUTMWGS84" Next

Selezionare la proiezione con id 2266 Next

Nuovo mapset :"nicola" Next

Riferirsi allo stesso Database "Datigrass" (e quindi cartella) creato in precedenza Next

Creare una nuova location "WGS84" Next

Selezionare la proiezione con id 2585 Next

Nuovo mapset :"nicola" Next

Aprire il primo mapset creato ToscanaGaussBoaga (coordinate Gauss-Boaga) Menù plugins -> Grass -> Apri mapset

menù Plugin / Grass / apri mapset Location: WGS84 Mapset: nicola

menù Plugin / Grass / Grass tools Vector / reproject vector ... / v.proj

I vettori sono riproiettati al volo nel sistema di coordinate del progetto.

Cambiare il sistema di coordinate del progetto da Gauss Boaga a WGS84

Se nel sistema WGS84 i due vettori non fossero allineati, controllare quale è il sistema di proiezione di base del layer in coordinate Gauss-Boaga. Legenda / tasto destro / proprietà / scheda generale / Cambia

Impostare il sistema di proiezione Gauss-Boaga

I vettori adesso dovrebbero risultare allineati

Ad ulteriore verifica che il file di proiezione prj non permette di definire correttamente gli scostamenti rispetto al sistema geografico WGS84, è possibile esportare come shp il vettore in coordinate Gauss-Boaga (tasto destro sulla legenda e "salva come shapefile ..."). Se si importa lo shape nella vista i punti risultano nuovamente sfalsati rispetto a quelli del vettore di origine.

A differenza del formato prj (ESRI), il più generco formato wkt, permette di specificare correttamente il rapporto tra un generico sistema di proiezione ed il sistema WGS84. Per creare un file wkt relativo ad una certa proiezione: Aprire un mapset di una location relativa alla proiezione di interesse (ad esempio ToscanaGaussBoaga / nicola) nella shell di grass digitare g.proj -p sono stampate sulla shell le caratteristiche della proiezione per stamparle in formato wkt digitare g.proj -w copiare la risposta in un file di testo vuoto Dare al file estensione wkt e lo stesso nome dello shape a cui si vuole associare

si imposta la regione di lavoro intorno all'alettrodotto appena digitalizzato con ampio margine W: 1674600 E: 1679700 S: 4847300 N: 4851900

si carica nella vista lo shapefile dell'edificato che sta nella cartella CARTOGRAFIA_BASE se si pianta tutto passare al punto successivo

Si tratta di uno shapefile molto pesante di cui a noi interessa solo una porzione.

dalla shell di GRASS lanciare il comando v.in.ogr dsn=C:/Archivio/Cartografia/SHP/CARTOGRAFIA_BASE layer=edificato_ab_2001_fi output=edificato spatial=1674600,4847300,1679700,4851900 l'opzione spatial=xmin,ymin,xmax,ymax permette di importare nel layer GRASS edificato soli gli elementi all'interno del rettangolo definito dalle coordinate di due suoi angoli Noi utilizziamo le stesse coordinate con cui abbiamo definito la regione di lavoro. E' necessario utilizzare l'opzione -o per riproiettare correttamente l'shp.

dal browser della toolbox GRASS si carica il layer 1 polygon dalla mappa GRASS edificato

la mappa dell'edificato è adesso molto più leggera e gestibile. Verificare che la tabella degli attributi abbia conservato i campi originari.

Creare un buffer di 400m intorno al tracciato della linea Toolbox GRASS / Modules / vector / Buffer / v.buffer

con la setessa procedura: Creare un buffer di 300m intorno al tracciato della linea

Impostare la regione di lavoro con l'estensione e la risoluzione desiderata per la griglia di punti

dalla shell di GRASS r.mapcalc surface=1

dalla shell di GRASS r.to.vect feature=point in=surface out=ptgriglia5m

Se la visualizzazione del grigliato appena generato fosse troppo lenta e più in generale per gestire la visualizzazione di mappe molto pesanti eliminare la spunta dalla check box "disegna". Ciò permette di effettuare varie operazioni (ad esempio la modifica della simbologia di un layer o dell'ordine di visualizzazione) senza che la vista sia ridisegnata. cambiare l'ordine di visualizzazione delle mappe aumentare lo zoom (la vista panoramica continua ad essere aggiornata) rimettere la croce nella check box "disegna" e visualizzare i punti della griglia

Toolbox GRASS / Modules / Vectors /extract features from vectors / v.select

input: ptgriglia5m input: buff300m output: ptgriglia5m_buff300m

eliminare la spunta dalla check box "disegna" cambiare l'ordine di visualizzazione delle mappe aumentare lo zoom (la vista panoramica continua ad essere aggiornata) impostare la trasparenza (50-60%) dei buffer rimettere la croce nella check box "disegna" e visualizzare i punti della griglia filtrati all'interno del buffer a 300m

rendere un layer grass puntuale a 3D in base al dtm da shell v.drape input=ptgriglia5m_buff300m type=point rast=dtm_090_10m method=bilinear output=ptgriglia5m3D_buff300m

esportare il layer grass come testo da shell v.out.ascii input=ptgriglia5m3D_buff300m format=point output=C:/Archivio/MasterGIS/Esercitazioni/Dati/grigliato_407_300m.txt fs="\t"

Questa parte del tutorial è solo esemplificativa-informativa e non si forniranno gli strumenti di modellazione e di calcolo. Vuole essere solo un esempio di come si procede con una delle applicazioni utilizzabili per la modellazione di elettrodotti ed il il calcolo dei campi.

esportare il layer sostegni_407_3d come file di testo formato point nella cartella DatiElettrodotti v.out.ascii input=sostegni_407_3d format=point output=C:/Archivio/MasterGIS/Esercitazioni/Dati/DatiElettrodotti/sostegni_407_3d.txt fs="\t"

In questo file sono riportate solo le coordinate dei sostegni ma non le altre informazioni (tipo e altezza) inserite al momento della digitalizzazione (tabella degli attributi).

Esporto la tabella dbf come file di testo echo "select * from sostegni_407_3d"|db.select fs=" " > C:/Archivio/MasterGIS/Esercitazioni/Dati/DatiElettrodotti/sostegni_407_3d.txt

Carico il file appena creato nel programma di modellazione

Assegno il tipo al primo sostegno (il portale)

Assegno il tipo al gruppo di sostegni 2-7. Completo l'assegnazione dei tipi come nei due casi precedenti, con riferimento alla tabella testuale.

Una volta assegnati i tipi, si modella il cablaggio della linea.

Importare il flie di testo appena creato nell'applicazione di calcolo del campo magnetico.

Calcolare il campo magnetico a diverse altezze da terra (quota del terreno definita dalla coordinata Z dei punti del grigliato)

Determinare la fascia di rispetto dei primo livello

Calcolare alcune curve isocampo a 1 m da terra considerando un modello del terreno semplificato

II risultati del calcolo sono salvati sia sotto forma di file di testo sia come shapefile. Caricare ad esempio lo shapefile relativo ad un calcolo su grigliato (ad esempio 1m da terra)

Modificare la simbologia dei punti nel modo indicato in figura ed aumentare lo zoom.

Caricare ad esempio lo shapefile relativo aL un calcolo su grigliato a 4m da terra. Modificare la simbologia dei punti nel modo indicato in figura.

Non si ha il diretto controllo sul colore associato al valore di campo.

Simile a "simbolo graduato" ma in questo caso la classe è costituita dai punti che assumono lo stesso valore del campo di classificazione. Modalità di visualizzazione adatta a casi in cui ci siano pochi valori del campo di classificazione che si ripetono per molti punti.

isolamento edificato con altezza superiore a 7 m: da shell v.extract input=edificato output=edificato7m where="heigth>7" --o isolamento punti di calcolo caratterizzati da campo maggore di 3 uT v.extract input=Bgriglia_7m output=Bgriglia_7m_3uT where="B>3" --o

Toolbox GRASS / Modules /Vector/Extract features from vector / v.select

estrazione dell'edificato on altezza superiore a 7 m con livelli di esposizione localmente maggiori di 3 microtesla in: edificato7m in: Bgriglia_7m_3uT out:edificato_7m_3uT

il modo più semplice per valutare il numero di abitanti esposti è quello di fare una query e copiare la risposta in un foglio di calcolo, in cui effettuare la somma echo "select ABITANTI01 from edificato_7m_3uT"|db.select

generazione di un raster che rappresenta i livelli di campo magnetico importare la griglia in GRASS Toolbox GRASS / Modules /import vector/v.in.ogr output: Bgriglia_7m

caricare il layer GRASS Bgriglia_7m utilizzare il comando GRASS Toolbox GRASS / Modules / raster /Interpolate surface / v.surf.idw per creare un raster che rappresenti i livelli di campo magnetico a 7 m dal terreno. Come già visto in precedenza è necessario utilizzare la shell di GRASS v.surf.idw input=Bgriglia_7m column="B" output=B_7m npoints=4

Toolbox GRASS / Modules /Import / Map Type conversion /v.to.rast Convert vector to raster using constant trasformare l'edificato a 7 m in raster

generare il nuovo raster in modo che se ci si trova dentro ad un edificio ed il livello di campo è maggiore di 3 uT allora la cella del raster vale 1. da shell: r.mapcalc 'ed_7m_3uT=if(edificato_7m==1&&B_7m>3,1,null())' r.mapcalc 'ed_7m_1to3uT=if(edificato_7m==1&&B_7m>1&&B_7m<3,1,null())' r.mapcalc 'ed_7m_04to1uT=if(edificato_7m==1&&B_7m>0.4&&B_7m<1,1,null())'

genero i raster in cui differenzio 4 regioni in base ai lielli di esposizione (copiare le seguenti righe e poi Maiusc+ clic tasto sinistro nella shell GRASS) r.mapcalc 'B_7m_3uT=if(B_7m>3,1,null())' r.mapcalc 'B_7m_1to3uT=if(B_7m>1&&B_7m<3,1,null())' r.mapcalc 'B_7m_04to1uT=if(B_7m>0.4&&B_7m<1,1,null())' r.mapcalc 'B_7m_0to0.4uT=if(B_7m>0&&B_7m<0.4,1,null())'

carico i 3 raster relativi ad esposizioni maggiori di 0,4 microtesla

trasformo i raster in vettori (poligoni) (copiare le seguenti righe e poi Maiusc+ clic tasto sinistro nella shell GRASS) r.to.vect input=B_7m_3uT output=B_7m_3uT feature=area --o r.to.vect input=B_7m_1to3uT output=B_7m_1to3uT feature=area --o r.to.vect input=B_7m_04to3uT output=B_7m_04to3uT feature=area --o visualizzare i vettori grass appena creati

Utilizzando il comando Overlay creo delle mappe vettoriali che rappresentano le porzioni di edificio esposte ad un certo intervallo di vlori di campo (http://grass.itc.it/grass63/manuals/html63_user/v.overlay.html) (copiare le seguenti righe e poi Maiusc+ clic tasto sinistro nella shell GRASS) v.overlay ainput=B_7m_3uT atype=area alayer=1 binput=edificato7m btype=area blayer=1 output=edificato_7m_3uT operator=and olayer=1,0,0 --o v.overlay ainput=B_7m_1to3uT atype=area alayer=1 binput=edificato7m btype=area blayer=1 output=edificato_7m_1to3uT operator=and olayer=1,0,0 --o v.overlay ainput=B_7m_04to1uT atype=area alayer=1 binput=edificato7m btype=area blayer=1 output=edificato_7m_04to1uT operator=and olayer=1,0,0 --o visualizzare i vettori grass appena creati utilizzando dei colori adatti ad evidenziare le elaborazioni fatte