Un breve richiamo: la Full Planar Topology
La Full Planar Topology è un modello topologico poligonale che soggiace ad alcune regole fondamentali:
- Connettività:
- gli archi sono tra loro connessi con un nodo in comune; ogni arco inizia e termina con un nodo (ogni arco ha due nodi); tutte le intersezioni tra archi sono nodi.
- Contiguità:
- poligoni adiacenti hanno un solo arco in comune; ogni arco divide due aree.
- Contenimento:
- un poligono chiuso da archi ha un'area misurabile; serie di archi definiscono perimetri e aree chiuse. Ogni area è circondata da archi e nodi. Ogni nodo è contenuto in aree.
Poste queste regole, il vettoriale polyTopology
della
location
simpleloc
rappresentato qui di seguito,
dovrebbe avere le tabelle riportate sotto:
![]() |
Le entità geometriche fondamentali: nodi, archi, poligoni. Il verso degli archi si riferisce al verso di digitazione. |
|
|
|
Obbiettivi
GRASS è un GIS topologico; le 3 tabelle fondamentali viste sopra devono quindi essere facilmente ottenibili.
Il comando
v.build map=polyTopology opt=dump
consente la visualizzazione di tutti gli elementi della topologia
del vettoriale
polyTopology
.
Lo scopo di questo esercizio è didattico: nella costruzione delle tre tabelle topologiche, si condurrà una simulazione senza ricorrere a v.build; questo consentirà di esplorare al meglio alcuni tools che interagiscono con le entità geometriche fondamentali (nodi,archi,poligoni).
Set di dati per l'esercitazione
La copertura poligonale scelta per l'esercitazione è polyTopology
,
contenuta nella location simpleloc.
Questa location è stata creata appositamente per fini didattici; per avere maggiori notizie vedere questo precedente blog.
La location simpleloc può essere scaricata da qui.
Procedura
Si creerà per prima la tabella Arc Topology, quindi la Node Topology e da ultimo la Polygon Topology.
Le fasi logiche sono le seguenti:
- creazione di un layer di archi (=boundaries) numerato univocamente;
- agli attributi degli archi vengono aggiunti gli identificativi (=cat) dei poligoni posti a destra e a sinistra di ciascun elemento;
- creazione di un layer di nodi (=punti) numerati univocamente;
- agli attributi dei nodi vengono aggiunti gli identificativi degli archi ai quali sono uniti;
- dalla combinazione delle tabelle dei nodi e degli archi si ottiene la Polygon Topology.
Arc Topology
I boundaries per default non hanno una categoria. Si crea un nuovo layer aggiungendo un valore di categoria univoco:
v.category in=polyTopology type=boundary \
layer=2 out=arcs option=add cat=1 --overwrite
Agli attributi del layer appena creato si aggiunge una serie di campi: alcuni conterrano solo testo e serviranno per la visualizzazione. Inoltre alcuni campi testuali offorno l'occasione per fare esercizio anche con l'SQL.
v.db.addtable arcs layer=2 columns="cat integer, \
arcLabel char(3), \
catPolyLeft integer, \
polyLeftLabel char(5), \
catPolyRight int, \
polyRightLabel char(5), \
arcLength double precision"
Si inizia a popolare il campo arcLabel
:
v.db.update map=arcs layer=2 column=arcLabel qcolumn="'a'||cat"
L'operatore SQL || concatena le stringhe; in questo caso concatena la stringa a con il contenuto del campo cat.
Con il tool v.to.db opt=sides si ottengono le categorie dei poligoni rispettivamente a sinistra e a destra di ciascun arco (boundary):
v.to.db arcs layer=2 option=sides col=catPolyLeft,catPolyRight
Si calcola inoltre la lunghezza di ciascun arco:
v.to.db arcs layer=2 option=length col=arcLength
Viene ora il momento di usare l'SQL: sulla base delle categorie delle aree ottenute con v.to.db opt=sides si popolano i due campi polyLeftLabel e polyRightLabel:
echo "UPDATE arcs_2 SET polyLeftLabel=IFNULL( \
(SELECT arcs.label FROM arcs \
WHERE arcs.cat=arcs_2.catPolyLeft),'universe')" | db.execute
echo "UPDATE arcs_2 SET polyRightLabel=IFNULL( \
(SELECT arcs.label FROM arcs \
WHERE arcs.cat=arcs_2.catPolyRight),'universe')" | db.execute
L'operazione effettuata merita qualche spiegazione; sul layer 1 del
vettoriale
arcs
vi sono gli attributi delle aree:
v.db.select arcs layer=1
cat|label|GRASSRGB
1|polyA|141:211:199
3|polyC|190:186:218
2|polyB|255:255:179
4|polyD|251:128:114
5|polyE|128:177:211
6|polyF|253:180:98
Nei campi catPolyLeft e catPolyRight del layer2 sono state memorizzate
le categorie delle aree
del layer 1: esiste quindi una relazione biunivoca fra i due campi del
layer 2 e il campo
cat del layer1. Grazie a questa relazione
si può stabilire un
join fra le due tabelle imponendo la
condizione WHERE
arcs.cat=arcs_2.catPolyLeft
e WHERE
arcs.cat=arcs_2.catPolyRight
.
Se tali condizioni non vengono soddisfatte significa che il valore
della cat del layer 1 non ha
corrispondenze nel layer 2: è il caso del poligono
universe (valore -1
), l'area
che, nel modello
topologico Full Planar Topology avvolge
idealmente ogni altro
oggetto.
Dopo queste operazioni gli attributi sul layer 2 sono:
v.db.select arcs layer=2
cat|arcLabel|catPolyLeft|polyLeftLabel|catPolyRight|polyRightLabel|arcLength
1|a1|3|polyC|1|polyA|6.242641
2|a2|2|polyB|1|polyA|3
3|a3|4|polyD|1|polyA|6.708204
4|a4|5|polyE|-1|universe|3.414214
5|a5|-1|universe|1|polyA|5.019765
6|a6|-1|universe|1|polyA|2.828427
7|a7|1|polyA|6|polyF|4.162278
8|a8|2|polyB|-1|universe|5.162278
9|a9|4|polyD|-1|universe|10.261297
10|a10|4|polyD|2|polyB|3
11|a11|5|polyE|-1|universe|5.828427
12|a12|5|polyE|6|polyF|2.414214
13|a13|-1|universe|6|polyF|7.071068
Si verifichi anche la situazione delle categorie sui due layers:
v.category arcs opt=report
Layer/table: 2/arcs_2
type count min max
point 0 0 0
line 0 0 0
boundary 13 1 13
centroid 0 0 0
area 0 0 0
all 13 1 13
Layer/table: 1/arcs
type count min max
point 0 0 0
line 0 0 0
boundary 0 0 0
centroid 6 1 6
area 0 0 0
all 6 1 6
A questo punto la costruzione della tabella Arc Topology è conclusa:
v.db.select arcs layer=2 \
col=arcLabel,polyRightLabel,polyLeftLabel
arcLabel|polyRightLabel|polyLeftLabel
a1|polyA|polyC
a2|polyA|polyB
a3|polyA|polyD
a4|universe|polyE
a5|polyA|universe
a6|polyA|universe
a7|polyF|polyA
a8|universe|polyB
a9|universe|polyD
a10|polyB|polyD
a11|universe|polyE
a12|polyF|polyE
a13|polyF|universe
Si procede ora alla visualizzazione :
./sl_draw.sh
d.font.freetype font=timesbd
d.vect -a arcs layer=1 display=shape
d.vect arcs layer=1 display=attr attrcol=Label llayer=1 \
color=black lsize=10 xref=center yref=bottom
d.vect arcs layer=2 display=shape,attr,dir attrcol=arcLabel \
type=boundary llayer=2 \
xref=center yref=center color=grey \
lsize=10 bgcolor=white bcolor=black
Lo script sl_draw.sh
può essere
scaricato da questo
indirizzo; deve essere ovviamente scomapttato, reso eseguibile (chmod
ugoa+x
sl_draw.sh
) e posto in una directory in modo che l'interprete
dei comandi possa
trovarlo (controllare con echo $PATH
).
![]() |
Il vettoriale arcs: sul layer 1 mantiene le aree; sul layer 2 gli archi (=boundaries) con gli attributi calcolati. |
where is sl_draw.sh script ?
RispondiElimina