Tutorial esercitazioni inquinamento elettromagnetico a bassa frequenza Aggiornato al 7 marzo 2009

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. Il tutorial è stato ulteriormente sviluppato nell'ambito delle edizioni 2007/2008 e 2008/2009 del master. Le esercitazioni seguono 2 ore di seminario sul tema dell'inquinamento elettromagnetico a bassa frequenza ed in particolare sulle linee elettriche ad alta tensione. In particolare nel 2009 è stata rilasciata la versione 1.0 di QGIS. 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 da meno tempo è possibile utilizzarlo anche in windows. Le versioni per windows, essendo meno mature sono anche anche meno stabili di quella per linux.
Si tratta in ogni caso di un programma in continua evoluzione e dalle grandi potenzialità e quindi vale la pena provare ad utilizzarlo. 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. Quanto detto per QGIS, vale a maggior ragione anche per GRASS, visto che il suo "porting" per windows è ancora più recente.
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 fino ad ora.

• Sito ufficiale QGIS
• Blog QGIS
• Download ultima versione QGIS
• Wiki QGIS: indice dei titoli
• Forum di QGIS
• Vecchio archivio mailing-list sviluppatori di QGIS (per controllare lo stato di sviluppo) (fino a marzo 2008)
• Vecchio archivio mailing-list utenti di QGIS(fino a marzo 2008)

• 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

Il file zip osgeopkgs.zip contiene qgis 1.0 e vari plugin per quest'ultimo, tra i quali grass. Utilizzando questo pacchetto è possibile installare qgis 1.0 anche su computer non connessi ad internet. Per fare questo è necessario: scompattare il pacchetto osgeopkgs.zip in c: usando l'opzione "use folder names" lanciare il file setup.bat

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

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 creare una nuova cartella C:/EsercitazioniQGIS2009/Qgis_proj dove salvare tutti i progetti di QGIS

menù Plugin / Gestione Plugins

selezionare GRASS e confermare (OK)

menù Plugin / Grass / New Mapset

Creare una nuova cartella "DatasetGrass" (ad esempio in C:/EsercitazioniQGIS2009/)

Selezionare la cartella appena creata e confermare

Next

Creare una nuova location "LocationToscanaGaussBoaga" Next

Selezionare la proiezione definita precedentemente Next

Impostare la regione GRASS sull'Italia Next

Nuovo mapset : Suggerimenti per l'assegnazione del nome: gli si potrebbe assegnare il proprio nome nel caso si utilizzi un Dataset/location accessibile da più utenti ("Nicola" ad esempio) nel caso di un pc a singolo utente potrebbe essere nominato in base alla particolare attività svolta ("MasterGIS" ad esempio). Next

Finish

menù Layer / aggiungi un raster

Cercare in C:\EsercitazioniQGIS2009\Dati\Raster e selezionare 275020.tif 275030.tif

menù Visualizza / vista massima

Se non si dovesse vedere nulla fare doppio clic sul nome del raster nella legenda Selezionare la scheda Generale

Cambiare il sistema di riferimento con quello generato precedentemente. NB: i raster sono immediatamente visualizzati se non si abilita la riproiezione al volo nelle proprietà del progetto.

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 Nella versione 0.9.2 e seguenti 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: Toolbox GRASS / vector / overlay /v.patch

Fare due prove, una selezionando ed una non selezionando la check box "Copy also attribute table" input=p275030,pc275030 output1=ppc275030noattributes output2=ppc275030

Nel secondo caso si ottiene una mappa puntuale a cui è associata una tabella degli attributi completa.

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.

Meglio se impostare le coordinate che definiscono i limiti della regione con dei multipli della risoluzione desiderata

Toolbox GRASS / Modules / raster / interpolate surface v.surf.idw chiamo il raster dtm_10m_12pt

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

Se il raster appare colorato uniformemente come nella precedente immagine fare doppio clic sul nome del raster nella legenda (oppure tasto destro proprietà) ed impostarle come nella seguente immagine (Stira a MinMax e poi Applica).

L'immagine sarà colorata in base alla quota attribuita a ciascun pixel secondo una mappa di colori in scala di grigi.

Ogni comando di GRASS, oltre che dalla toolbox, può essere lanciato utilizzando la shell di grass cioé una finestra di comandi. Si crea un nuovo dtm utilizzando un numero diverso di punti (4 in questo caso) nella GRASS shell v.surf.idw input=ppc275030 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_idw12pt_idw4pt=dtm_10m_12p-dtm_10m_4p'

Se si vogliono visualizzare correttamente i raster appena generati si deve modificarne le propreià, come fatto per il primo (doppio clic sulla legenda-> proprietà ->Scheda Simbologia -> Stira MinMax -> Apply)

Visualizzare l'istogramma di uno dei dtm, ad esempio dtm_10m_12pt. Tasto destro / proprietà / istogramma / tasto aggiorna. Si nota che le quote sono comprese tra 28 e 153 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_12pt 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_4pt color=rules

E' possibile lanciare il comando dall'interfaccia grafica del gestore strumenti di GRASS: Toolbox Grass -> raster -> r.colors.rules in tal caso eliminare end finale dal file delle regole

E' necessario eliminare end finale dal file delle regole altrimenti si ottiene il seguente errore

Dopo aver eliminato end finale dal file delle regole altrimenti il comando va a buon fine

E' inoltre possibile utilizzare il comando r.colors utilizzando "valori assoluti" r.colors map=dtm_10m_12pt 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.

Si esegue l'analisi per le pendenze (r.slope.aspect.slope)

Si esegue l'analisi per le pendenze (r.slope.aspect.aspect)

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").

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 sost407 apri strumenti GRASS / v.in.ogr eliminare il layer sostegni_407.shp caricare il layer grass sost407 (icona o menù Plugin / GRASS / Aggiungi vettoriale GRASS)

Aggiunta etichette ai sostegni evidenziare sost407 nella legenda tasto destro / proprietà oppure doppio clic sul nome del layer Scheda etichette: sub-scheda posizione: campo con etichetta: codice spostamento unità: punti

La versione 1.0.1 di QGIS ha un baco che non permette di modificare lo stile delle etichette Nella seguente immagine si riportano i risultati ottenuti con le versioni precedenti.

caricare il raster GRASS dtm_10m_12pt creato in precedenza

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

input: sost407 output: sost407_3d raster: dtm_10m_12pt 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

si usa la shell v.to.db map=sost407_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=sost407_3d format=point output=C:/EsercitazioniQGIS2009/WORK/s407_point.txt fs="\t"

aprire con wordpad il file di testo appena creato 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:/EsercitazioniQGIS2009/WORK/s407_std.txt output=tracciato407 --overwrite format=standard

caricare il layer vettoriale tracciato407

NOTA BENE APPROCCIO NON FUNZIONANTE CON QGIS 1.0.1: LA SHELL DI GRASS NON RISPONDE COME DOVREBBE NEL CASO IL NUMERO DEGLI ARGOMENTI NON SIA QUELLO DESIDERATO ED INOLTRE QUANDO SI APRE IL LAYER GRASS GENERATO QGIS VA IN CRASH E AL SUCCESSIVO RIAVVIO NON SI RIESCE AD APRIRE IL MAPSET per sbloccare il mapset eliminare il file .gislock nella radice della cartella del mapset 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 v.to.lines ed incollarlo nella cartella C:\Programmi\OSGEO\apps\grass\grass-6.4.0RC3\scripts Copiare il file v.to.lines.bat ed incollarlo nella cartella C:\Programmi\OSGEO\apps\grass\grass-6.4.0RC3\bin Lo script richiede che si specifichi: il nome del layer puntuale di ingresso il nome del layer delle linee in uscita Una volta copiato lo script provare a digitare nella shell: v.to.lines input=sost407_3d output=tracciato407bis

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 C:\Programmi\OSGEO\apps\qgis\grass\modules modificare il file C:\Programmi\OSGEO\apps\qgis\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 ED50/utm ZONE 32N Next

Nuovo mapset :"MasterGIS" 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 WGS84/UTM zone 32N Next

Nuovo mapset :"MasterGIS" Next

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

Creare una nuova location "WGS84" Next

Selezionare la proiezione WGS84 Next

Nuovo mapset :"MasterGIS" Next

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

importare con v.in.ogr lo shapefile C:\EsercitazioniQGIS2009\Dati\DatiElettrodotti\shp\sostegni_40700.shp

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

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

I vettori sono riproiettati al volo nel sistema di coordinate del progetto (sempre che sia attivata l'opzione della riproiezione al volo nelle proprietà 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.

Toolbox GRASS -> v.in.ogr

selezionare l'opzione "Limit import to current region"

oppure dalla shell di GRASS lanciare il comando v.in.ogr dsn=C:\Archivio\MasterGIS\Esercitazioni\Dati\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 Con le vecchie versioni del programma Toolbox GRASS / Modules / vector / Buffer / v.buffer

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

Con la nuova versione del programma è necessario usare la shell: v.buffer input=tracciato407 output=buffer400 distance=400 v.buffer input=tracciato407 output=buffer300 distance=300

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

oppure con la Toolbox GRASS

esportare il layer grass come testo da shell v.out.ascii input=ptgriglia5m3D_buff300m format=point output=C:\EsercitazioniQGIS2009\WORK\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:\EsercitazioniQGIS2009\WORK\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

OPPURE Toolbox GRASS / Modules /Vector/Spatial analysis /v.extract.where

importazione grigliato a 7 m

isolamento punti di calcolo caratterizzati da campo maggore di 3 uT da shell v.extract input=Bgriglia_7m output=Bgriglia_7m_3uT where="B>3" --o

OPPURE, come in precedenza, Toolbox GRASS / Modules /Vector/Spatial analysis /v.extract.where

ESTRAZIONE estrazione dell'edificato Con altezza superiore a 7 m con livelli di esposizione localmente maggiori di 3 microtesla: 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 se non lo si è già fatto 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. Dalla shell di GRASS v.surf.idw input=Bgriglia_7m column="B" output=B_7m npoints=4

Oppure dalla Toolbox 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())"

modifico le mappe di colori associate ai 3 raster da shell: r.colors map=ed_7m_3uT color=rules 1 255 0 0 end r.colors map=ed_7m_1to3uT color=rules 1 0 0 255 end r.colors map=ed_7m_04to1uT color=rules 1 0 255 0 end

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_04to1uT output=B_7m_04to1uT 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

Importare come layer GRASS, con v.in.ogr i seguenti shp: C:\EsercitazioniQGIS2009\Dati\DatiElettrodotti\shp\sostegni_40700_strutturasostegni.shp C:\EsercitazioniQGIS2009\Dati\DatiElettrodotti\shp\catenarie_40700_catenarie.shp C:\EsercitazioniQGIS2009\Dati\DatiElettrodotti\shp\sostegni_40900_strutturasostegni.shp C:\EsercitazioniQGIS2009\Dati\DatiElettrodotti\shp\shp chiamandoli rispettivamente: strutsost407 cat407 strutsost409 cat409 ed avendo cura di importarli come layer 3D (si veda la seguente immagine

Nella ToolboxGrass, aprire il modulo NVIZ

Utilizzare come dtm uno dei raster creati durante le esercitazioni e scegliere uno dei vettori grass 3D appena importati

L'estensione della superficie visualizzata è quella dell'area di lavoro grass attualmente impostata. Con la freccia verde è possibile modificare posizione ed orientamento del punto di vista.

Menù visualize -> Vector lines / 3D polygons New aggiungere i layer delle catenarie e delle strutture dei sostegni di entrambe le linee

Per visualizzare l'edificato è prima necessario "estruderlo" con il comando v.extrude v.extrude input=edificato output=edificato_vol elevation=dtm_10m_12pt hcolumn=heigth

Torno in NVIZ Menù visualize -> Vector lines / 3D polygons New aggiungere il layer edificato_vol