Come risolvere gli allineamenti per intersezione

Gianni Rossi

Gli allineamenti per intersezione sono uno artificio gettonatissimo nei rilievi topografici e catastali eseguiti con tecnologia GPS perché risolvono egregiamente (e comodamente) il problema del mancato ricevimento del segnale satellitare nella rilevazione degli spigoli di fabbricato quando l’edificio stesso (o un agglomerato di più edifici) ne impedisce la ricezione. Lo schema è molto semplice, si rilevano con il GPS due punti nei pressi dello spigolo (posti ad una distanza tale da ricevere il segnale) e da questi si misura, con il disto o la cordella metrica, la distanza (orizzontale) allo spigolo. L’intersezione (da cui il nome) delle due circonferenze determina la posizione dello spigolo da rilevare.

Operativamente, dunque, è tutto estremamente semplice. La risoluzione topografica è invece un po’ più complicata, cercherò qui di illustrarla a beneficio dei colleghi più appassionati che non disdegnano di capire anche il “come” le misure prese in campagna producono i risultati attesi. Lo schema geometrico è il seguente:

Figura 1Lo schema dell’intersezione tra due circonferenze. Legenda:
C1(a, b) e C2(c, d) = centri e rispettive coordinate
D = distanza tra i due centri
P1 , P2 = punti di intersezione
A = area del triangolo formato dai due centri e i punti di intersezione: C1-C2-P1 o C1-C2-P2

Per risolverlo bisogna per prima cosa verificare se effettivamente le due circonferenze si intersecano. Questo avviene solo se la distanza D tra i due centri è maggiore della differenza tra i due raggi e minore della loro somma. Deve cioè essere rispettata la condizione:

\[r_1-r_2 ≤ D ≤ r_1+r_2 \]

Perché questo? Lo vediamo esaminando le cinque possibili geometrie:

  1. Se la differenza tra i due raggi è maggiore di D, la seconda circonferenza risulta completamente inglobata alla prima.
  1. Se la differenza tra i due raggi è esattamente uguale a D, le due circonferenze sono tangenti tra loro, con la seconda completamente interna alla prima.
  1. Se la somma tra i due raggi è esattamente uguale a D, le due circonferenze sono tangenti tra loro, con la seconda completamente esterna alla prima.
  1. Se la differenza tra i raggi è minore di D e la somma è invece maggiore, allora si intersecano in due punti.
  1. Se la somma tra i due raggi è minore di D, le due circonferenze sono completamente separate ed esterne una dall’altra.

La condizione 4 è pertanto quella da verificare per capire se sussistono gli estremi per la risoluzione. Se così, il calcolo inizia con il trovare l’area A di uno dei due triangoli (uguali tra loro) formati dai due centri e il rispettivo punto di intersezione tra le circonferenze, vale a dire il triangolo C1-C2-P1 oppure C1-C2-P2. Questa area si trova con la seguente formula[1]:

\[A=\frac{\sqrt{(D+r_1+r_2 )(D+r_1-r_2 )(D-r_1+r_2 )(-D+r_1+r_2 )}}{4} \]

Dopodiché, con questo dato, più la distanza D tra i due centri e le coordinate degli stessi (a, b, c, d), le coordinate dei due punti di intersezione si trovano con queste altre formule:

\[ x_{P1,P2}=\frac{a+c}{2}+\frac{(c-a)(r_1^2-r_2^2)}{2 D^2} ± 2\frac{b-d}{D^2} A \]
\[ y_{P1,P2}=\frac{b+d}{2}+\frac{(d-b)(r_1^2-r_2^2)}{2 D^2} ∓ 2\frac{a-c}{D^2} A \]

I simboli ± e che precedono la parte a destra delle formule sono da intendersi così:
Per le x:
applicando il + si trova la coordinata di P1;
applicando il – si trova la coordinata di P2.
Per le y:
applicando il – si trova la coordinata di P1;
applicando il + si trova la coordinata di P2.

Naturalmente si deve poi stabilire quale dei due punti è quello desiderato. Questo lo si desume dalla convenzione standard usata in topografia, e cioè: mi posiziono in uno dei due punti rilevati (i centri C1 o C2 delle circonferenze), guardo l’altro punto rilevato e vedo se il punto determinato dalle due distanze è alla mia sinistra (devo girarmi in senso antiorario) oppure alla mia destra (devo girarmi in senso orario). Ad esempio, nel caso qui esposto mi posiziono in C1 e guardo C2, quindi:

  • Il punto alla mia sinistra è P1;
  • Il punto alla mia destra è P2.

Questa convenzione è stata fatta propria anche da Pregeo (vedi l’articolo Come rilevare per allineamenti punti non visibili o non accessibili) che la codifica mediante l’inserimento dell’angolo fittizio posto nelle due righe dei due allineamenti che definiscono l’intersezione. A questo angolo si può dare un valore qualsiasi, anche se si usa inserire 50g, cioè metà di un angolo retto, perché questo è il valore al quale dovrebbe avvicinarsi l’angolo effettivo del triangolo formato in campagna[2]. L’importante è il segno che si antepone al valore:

  • +50g se il punto è alla destra (orario) rispetto alla direzione origine-orientamento dell’allineamento;
  • -50g se il punto è alla sinistra (antiorario) rispetto alla direzione origine-orientamento dell’allineamento.

Lo vedremo sviluppando l’esempio numerico del rilievo di Figura 2 qui sotto. Si tratta di un rilievo GPS (tabella gialla al centro) che ha determinato i due punti 1001 e 1002 dai quali si è poi definito il punto 1 (spigolo di fabbricato non stazionabile) mediante le rispettive distanze riportate nella colonna Dist. della tabella verde degli allineamenti (in alto):

Figura 2Esempio di intersezione: dai punti GPS 1001 e 1002 (tabella gialla al centro) sono definite le due circonferenze di raggio indicato nella colonna “Dist.” della tabella verde degli allineamenti in alto. La tabella in basso riporta i risultati del calcolo.

Dalla tabella del calcolo (in basso in Figura 2), riscontriamo le coordinate dei due punti 1001 e 1002 derivanti dal rilievo GPS. Ignoriamo invece per il momento le coordinate del punto 1 che sono proprio lo scopo di questo esercizio. La prima operazione da fare è il calcolo della distanza D tra i punti 1001 e 1002, vale a dire i centri delle due circonferenze. Questo risultato si ottiene semplicemente applicando il teorema di Pitagora alla differenza di coordinate dei punti stessi:

\[ D=\sqrt{([216.823-210.199]^2+[-71.231-(-76.880)]^2}=8.706 \]

Verifichiamo ora la condizione di intersezione delle circonferenze data dal confronto tra la distanza  appena trovata e la differenza / somma dei due raggi:

\[ 7.326-5.963 < 8.706 < 7.326+5.963 \]
\[ 1.363 < 8.706 < 13.289 \]

Soddisfatta questa verifica, procediamo con il calcolo dell’area A dei due triangoli formati dai raggi e da ciascuno dei due punti di intersezione, sviluppando la formula iniziale vista sopra. Troviamo dapprima i quattro fattori della moltiplicazione che appare sotto radice quadrata al numeratore della frazione:

\[ D+r_1+r_2= 8.706 +5.963+7.236=21.905 \]
\[ D+r_1-r_2 = 8.706 +5.963-7.236=7.433 \]
\[ D-r_1+r_2 = 8.706-5.963+7.236=9.979 \]
\[ -D+r_1+r_2= -8.706 +5.963+7.236=4.493\]

Dopodiché troviamo l’area A:

\[ A=\frac{\sqrt{21.905∙7.433∙9.979∙4.493}}{4}=21.360 m^2 \]

Con questo dato, più le coordinate dei due centri 1001 e 1002, siamo ora in grado di calcolare le coordinate dei due punti di intersezione con le formule viste sopra. Troviamo dapprima i tre addendi di ciascuna formula e poi le coordinate dei due punti di intersezione.
Per le x:

\[ \frac{a+c}{2}=\frac{210.199+216.823}{2}=213.511 \]
\[ \frac{(c-a)(r_1^2-r_2^2 )}{2 D^2}=\frac{(216.823-210.199)∙(5.963^2-7.236^2)}{2∙ 8.706^2}=-0.734 \]
\[ 2\frac{b-d}{D^2} A = 2\frac{-76.880-(-71.231)}{8.706^2} 21.360 = -3.184 \]
\[ x_{P1}=213.511+(-0.734) + (-3.184)=209.593 \]
\[ x_{P2}=213.511+(-0.734)- (-3.184)=215.960 \]

Per le y:

\[ \frac{b+d}{2}=\frac{-76.880+(-71.231)}{2}=-74.055 \]
\[ \frac{(d-b)(r_1^2-r_2^2 )}{2 D^2}=\frac{[-71.231-(-76.880)]∙(5.963^2-7.236^2)}{(2∙8.706^2}=-0.626 \]
\[ 2\frac{a-c}{D^2}A=2 \frac{210.199-216.823}{8.706^2} 21.360=-3.734 \]
\[ y_{P1}=-74.055+(-0.626)- (-3.734)=-70.948 \]
\[ y_{P2}=-74.055+(-0.626) + (-3.734)=-78.416 \]

I due punti di intersezione hanno pertanto queste coordinate:

\[ P_1=209.593,-70.948 \]
\[ P_2=215.960,-78.416 \]

Dobbiamo ora stabilire quale dei due è il punto cercato. È evidente che se disegniamo sul CAD sia i punti 1001 e 1002 che i punti di intersezione P1 e P2 appena calcolati, possiamo appurare immediatamente quale dei due è il punto corretto: basterà verificare quello che rispetta le direzioni oraria e antioraria definita per i due allineamenti 1001-1002 e 1002-1001 dagli angoli fittizi +50 e -50. Ma questo è ciò che possiamo fare noi manualmente, mentre invece un algoritmo deve risolvere analiticamente l’ambiguità.

Come si fa?

Beh, esiste più di un procedimento, quello che descrivo qui è la procedura adottata in Geocat. Con riferimento alla Figura 3, consiste nel trovare gli “angoli interni con segno” dei due triangoli, cioè i valori angolari muniti di segno positivo o negativo a seconda che, per portare la congiungente 1001-1002 sopra la congiungente tra ciascuno di questi due punti e il punto di intersezione, si debba girare in senso orario oppure in senso antiorario. Trovati questi “angoli interni con segno”, è poi sufficiente confrontare in quale dei due triangoli i segni sono uguali a quelli degli angoli fittizi +50 e -50 degli allineamenti. Vediamo i vari passaggi da compiere. Gli angoli interni (con o senza segno) dei triangoli in questione si trovano per differenza degli gli azimut tra il vertice considerato e gli altri due vertici del triangolo.

Figura 3Una volta calcolati i due punti di intersezione, bisogna stabilire quale dei due è quello cercato. Per fare questo, si può adottare il calcolo degli “angoli interni con segno” dei due triangoli 1001-P1-1002 e 1001-P2-1002.

Per prima cosa troviamo l’azimut tra i punti 1001 e 1002 e gli azimut tra ciascuno di questi due punti e i due punti di intersezione. L’azimut ϑPQ tra due punti P e Q di coordinate note si trova calcolando dapprima l’arcotangente del rapporto tra la differenza delle x e quella delle y. Il valore assoluto così trovato corrisponde all’angolo al vertice e va poi rapportato al valore effettivo in funzione del quadrante in cui si trova il punto Q rispetto al punto P. Quest’ultima operazione è illustrata nella Figura 7 dell’articolo Come risolvere l’ambiguità della stazione libera.

\[ ϑ_{P-Q}=arctan\frac{x_q-x_p}{y_q-y_p} →rapporto\hspace{2mm}al\hspace{2mm}quadrante \]

Si deve però fare attenzione al fatto che, nel calcolare la tangente di triangoli rettangoli con lati così corti (pochi metri), è necessario utilizzare molti decimali, altrimenti l’angolo calcolato risulta troppo approssimato. Nei calcoli che seguono le lunghezze sono riportate con soli 3 decimali ma in realtà ne sono stati utilizzati ben 6. I risultati possono pertanto essere leggermente diversi da quelli ottenuti riproducendoli con la calcolatrice. Questa differenza non ha tuttavia alcuna rilevanza perché il valore degli angoli ci interessa unicamente per stabilirne il segno. Procediamo a calcolare tutti gli azimut necessari, iniziando dal 1001-1002:

\[ ϑ_{1001-1002}=arctan\frac{216.823-210.199}{-71.231-(-76.880)}=55.0520 \]
\[ ϑ_{1001-1002}→1°Quadrante=55.0520 \]

L’azimut reciproco 1002-1001 si calcola, come noto, aggiungendo 200g se l’azimut 1001-1002 è inferiore a 200g, oppure sottraendo 200g se è superiore:

\[ ϑ_{1002-1001}=55.0520+200=255.0520 \]

Passiamo ora agli azimut ϑ1001-P1, ϑ1001-P2, ϑ1002-P1, ϑ1002-P2:

\[ ϑ_{1001-P1}=arctan\frac{209.594-210.199}{-70.948-(-76.880)}=6.4803 \]
\[ ϑ_{1001-P1}→4°Quadrante = 400-6.4803=393.5197 \]
\[ ϑ_{1001-P2}=arctan\frac{215.960-210.199}{-78.416-(-76.880)}=83.4156 \]
\[ ϑ_{1001-P2}→2°Quadrante = 200-83.4156=116.5844 \]
\[ ϑ_{1002-P1}=arctan\frac{209.593-216.823}{-70.948-(-71.231)}=97.5040 \]
\[ ϑ_{1002-P1}→4°Quadrante = 400-97.5040=302.4960 \]
\[ ϑ_{1002-P2}=arctan\frac{215.960-216.823}{-78.412-(-71.231)}=7.6080 \]
\[ ϑ_{1002-P2}→3°Quadrante = 200+7.6080=207.6080 \]

Con gli azimut appena trovati passiamo ora a calcolare gli “angoli interni con segno”. Questo calcolo si esegue facendo semplicemente la differenza tra i due azimut interessati. Se questo valore è compreso nell’intervallo tra -200g e +200g, abbiamo già il valore cercato. Viceversa, se fuoriesce da tale range, va rapportato a 400g, come segue:

\[ α=ϑ_1-ϑ_2 \]
\[ se→ α>+200→α=α-400 \]
\[ se→ α<-200→α=α+400 \]

Procediamo:

\[ ang\hspace{2mm}P_1-1001-1002=ϑ_{1001-P1}-ϑ_{1001-1002} \]
\[ ang\hspace{2mm}P_1-1001-1002=393.5197-55.0520=338.4677 \]
\[ ang\hspace{2mm}P_1-1001-1002\hspace{2mm}→\hspace{2mm}>+200\hspace{2mm}→\hspace{2mm}338.4677-400=-61.5323 \]
\[ ang\hspace{2mm}P_2-1001-1002=ϑ_{1001-P2}-ϑ_{1001-1002} \]
\[ ang\hspace{2mm}P_2-1001-1002=116.5844-55.0520=61.5323 \]
\[ ang\hspace{2mm}P_2-1002-1001=ϑ_{1002-P1}-ϑ_{1001-1002} \]
\[ ang\hspace{2mm}P_2-1002-1001=302.4960-255.0520=47.4440 \]
\[ ang\hspace{2mm}P_2-1002-1001=ϑ_{1002-P2}-ϑ_{1002-1001} \]
\[ ang\hspace{2mm}P_2-1002-1001=207.6080-255.0520=-47.4440 \]

Analizziamo i risultati degli angoli interni in 1001 e 1002 appena calcolati e verifichiamo che i valori corrispondenti agli angoli fittizi di +50 e -50 con i quali avevamo definito la direzione del punto di intersezione sono quelli del triangolo P2-1002-1001, ovvero con positivo l’angolo in 1001 (61.5323) e negativo l’angolo in 1002 (-47.4440). Il punto di intersezione cercato è pertanto P2. La conclusione dell’esercizio sancisce che le coordinate del punto 1 misurato per doppia circonferenza dai punti GPS 1001 e 1002 sono:

\[ x_1=215.960\hspace{22mm}y_1=-78.416 \]

Le stesse coordinate sono infatti quelle risultanti dal calcolo di Geocat mostrato in Figura 2 (tabella in basso).

Bene, abbiamo trovato le coordinate X-Y del punto cercato risolvendo (analiticamente) anche l’ambiguità del doppio punto di intersezione.

E la quota del punto?

In effetti l’algoritmo che abbiamo visto è puramente planimetrico perché misura due distanze orizzontali senza alcuna informazione sul dislivello tra i due punti generatori dei quali invece, essendo stati rilevati da GPS o TS, è nota la quota. Naturalmente se si trattasse di un rilievo topografico in cui l’altimetria fosse un requisito richiesto, questo artificio non sarebbe assolutamente sufficiente. Nei rilievi in ambito catastale o per la ricostruzione di confini la precisione altimetrica non è in genere richiesta. Tant’e che Pregeo stesso permette di attribuire ai punti determinati per allineamenti mediante le pseudo-livellazioni da un estremo o dal mezzo. Vengono definite “pseudo-livellazioni” perché non vengono svolte con tutti i crismi e la precisione delle livellazioni vere e proprie. Anzi, normalmente si tratta di una semplice stima a vista fatta dal tecnico del dislivello tra i punti effettivamente rilevati (generatori dell’allineamento) e il punto determinato per intersezione. Tuttavia, se si desidera comunque ottenere una buona precisione anche sulla quota del punto si possono applicare le pseudo-livellazioni da un estremo o dal mezzo in maniera accurata. Per quella dal mezzo, ad esempio, con riferimento alla Figura 4 si opera come segue:

  • si sistema una palina, o l’asta del prisma, sul punto di cui si conosce la quota (101 in Figura 4, rilevato da GPS);
  • da questo, con un disto o con la cordella metrica, si prende la misura orizzontale allo spigolo e quella inclinata sul punto a terra (102 in Figura 4);
  • date le due distanze si calcola il dislivello con il teorema di Pitagora;
  • si misura l’altezza da terra sulla palina;
  • si eseguono i calcoli riportati in Figura 4 a seconda che il punto si trovi in alto o in basso.

Figura 4Una delle possibilità di sfruttare la pseudo-livellazione da un estremo per determinare la quota di uno spigolo di fabbricato determinato indirettamente (per allineamenti) a partire da un punto rilevato con il GPS nelle sue vicinanze. A sinistra l’esempio con lo spigolo in basso, a destra con lo spigolo in alto.

Se invece si adotta la pseudo-livellazione dal mezzo, lo schema è quello di Figura 5 che mostra come sfruttare questa tecnica per attribuire la quota ad uno spigolo di fabbricato (102) a partire dal punto nelle sue vicinanze rilevato con il GPS (101). I passaggi sono del tutto analoghi a quelli già visti sopra per la livellazione da un estremo e hanno poco a che fare con una livellazione. In pratica si tratta semplicemente di stimare il dislivello tra i due punti.

Figura 5Nella pseudo-livellazione dal mezzo il livello fittizio viene immaginato tra il punto noto (rilevato da TS o GPS, 101) e quello al quale si vuole attribuire la quota (spigolo, 102).


[1] Sarebbe anche interessante dimostrare sia questa formula che le successive per il calcolo delle coordinate, ma lo sarebbe per gli appassionati di matematica applicata alla topografia (come il sottoscritto), più che per gli appassionati di sola topografia. Si tratta infatti di formule che applicano il calcolo vettoriale la cui dimostrazione risulterebbe pertanto incomprensibile a chi non conosce questa branca della matematica-geometria. Tuttavia, su uno dei forum nei quali ho segnalato questo articolo (Facebook e altri social) un appassionato geometra topografo mi ha chiesto ripetutamente da dove salta fuori questa formula:

\[A=\frac{\sqrt{(D+r_1+r_2 )(D+r_1-r_2 )(D-r_1+r_2 )(-D+r_1+r_2 )}}{4} \tag{$1$} \]

Ho quindi pensato di fargli cosa gradita nel darne qui di seguito una dimostrazione che, pur evitando di entrare nel calcolo vettoriale, possa comunque riportarla nell’ambito della geometria studiata a scuola di geometra. La formula può infatti essere ricondotta a quella molto più nota di Erone:

\[A=\sqrt{P(P-r_1)(P-r_2)(P-D)} \tag{$2$} \]

Dove P è il semi-perimetro:

\[P=\frac{r1+r2+D}{2} \]

Vediamo perché le due formule (1) e (2) sono equivalenti analizzando i quattro fattori al numeratore (sotto radice quadrata) della formula (1) Il 1o fattore non è che il perimetro del triangolo:

\[D+r_1+r_2=2P \]

e quindi:

\[\frac{D+r_1+r_2}{2}=P \tag{$3$} \]

Il 2o fattore si può riscrivere così:

\[D+r_1-r_2=(D+r_1+r_2)-2r_2 \]

Dove gli addendi tra parentesi sono nuovamente il perimetro 2P del triangolo, e dunque si ha:

\[D+r_1-r_2=2P-2r_2 \]

e dividendo ambo i membri per 2:

\[\frac{D+r_1-r_2}{2}=P-r_2 \tag{$4$}\]

Analogamente per il 3o e 4o fattore si ha:

\[\frac{D-r_1+r_2}{2}=P-r_1 \tag{$5$}\]
\[\frac{-D+r_1+r_2}{2}=P-D \tag{$6$}\]

Sostituiamo ora le quattro espressioni (3), (4), (5), (6) nella formula di Erone (2):

\[A=\sqrt{\frac{D+r_1+r_2}{2}\hspace{2mm}\frac{D+r_1-r_2}{2}\hspace{2mm}\frac{D-r_1+r_2}{2}\hspace{2mm}\frac{-D+r_1+r_2}{2}} \]

Moltiplichiamo le frazioni sotto radice:

\[A=\sqrt{\frac{(D+r_1+r_2 )(D+r_1-r_2 )(D-r_1+r_2 )(-D+r_1+r_2 )}{16}} \]

Eliminiamo la radice quadrata al denominatore (calcolandola) e otteniamo la formula (1):

\[A=\frac{\sqrt{(D+r_1+r_2 )(D+r_1-r_2 )(D-r_1+r_2 )(-D+r_1+r_2 )}}{4} \]

Sono quasi certo che, a questo punto, il collega che mi ha chiesto di spiegare questa formula mi dirà: Ma allora perché non hai usato la formula di Erone? La domanda è logica, ma lo è anche la risposta: il motivo è che l’intero algoritmo applica il calcolo vettoriale e come tale ho applicato anche all’area tale approccio.

[2] Infatti, per gli allineamenti con i quali si definiscono i vertici di un fabbricato ortogonale a questo angolo viene dato il valore fittizio di 100g.


%d blogger hanno fatto clic su Mi Piace per questo: