Sensori Passivi e Attivi per la Rappresentazione 3D delle Città PDF

CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città SENSORI PASSIVI E ATTIVI AVIOTRASPORTATI PER LA RAPPRESENTAZIONE 3D DELLE CITTÀ 4.1. INTRODUZIONE Negli ultimi anni si è diffuso sempre più l'uso di sensori aerei di tipo attivo e/o passivo per applicaz...

CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città SENSORI PASSIVI E ATTIVI AVIOTRASPORTATI PER LA RAPPRESENTAZIONE 3D DELLE CITTÀ 4.1. INTRODUZIONE Negli ultimi anni si è diffuso sempre più l'uso di sensori aerei di tipo attivo e/o passivo per applicazioni di fotogrammetria e telerilevamento (Photogrammetry & Remote Sensing - PRS). Questo successo è dovuto principalmente alla flessibilità e alla capacità di raccogliere dati geometrici ad alta risoluzione, anche grazie all'impressionante diffusione di Unmanned Aerial Systems (UAS) nel campo della PRS. L'utilizzo di tali piattaforme versatili e flessibili ha portato ad una nuova e importante rivitalizzazione dell'acquisizione di geo-dati da piattaforme aeree. Nonostante sia importante sottolineare che i rilievi con e senza equipaggio presentano differenze in diversi aspetti quali la durata del volo, la copertura a terra e le tecniche di acquisizione dei dati, condividono un background comune. Uno degli aspetti principali è senza dubbio l'uso di telecamere digitali per l'acquisizione di immagini, anche se le obsolete telecamere analogiche sugli aerei con equipaggio possono essere ancora impiegate in progetti minori. È consuetudine dividere i sensori in passivi e attivi (Gomarasca, 2009). Come noto, i sensori ottici passivi possono rilevare l'energia elettromagnetica naturale (radiazione) che viene riflessa dall'oggetto osservato. Infatti, a seconda dei requisiti dell'applicazione desiderata, esistono sensori con diverse geometrie di acquisizione dati, formato, risoluzione geometrica e radiometrica/spettrale (compreso il numero, la lunghezza d'onda centrale e l'ampiezza delle bande spettrali). I sensori attivi si basano sull'emissione di segnali elettromagnetici che vengono riflessi dalla superficie da ricostruire o analizzare. Il sensore stesso può registrare il ritorno del segnale da elaborare. I principali tipi di tecnologie di sensori attivi adatti ad essere installati su piattaforme aeree sono i seguenti: (i) LiDAR (Light Detection and RAnging o Airborne Laser Scanning (ALS), e (ii) sensori radar a microonde. Va anche menzionato che altre tecnologie sono state sperimentate sperimentalmente su piattaforme UAS (per esempio, le fotocamere Time-of-Flight, mentre nuovi sensori come i LiDAR a singolo fotone sono ora in rapido sviluppo e si prevede che possono acquisire una grande quota di mercato nel prossimo futuro. CORSO DI GIS E 3D CITY MODEL 60  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città 4.2. SENSORI PASSIVI AVIOTRASPORTATI: LE FRAME CAMERA 4.2.1. CONCETTI DI BASE La fotogrammetria è quella tecnica che consente di definire la posizione, la forma e le dimensioni degli oggetti sul terreno, utilizzando le informazioni contenute in opportune immagini fotografiche degli stessi oggetti, riprese da punti diversi. Per poter determinare le posizioni dei punti di un oggetto nell’ambiente reale del territorio utilizzando le posizioni dei punti corrispondenti sulla fotografia, è necessario definire le relazioni geometriche fra le posizioni tridimensionali dei punti dell’oggetto e quelle delle loro immagini sul piano della fotografia. A ogni punto dell’oggetto tridimensionale (spazio oggetto) corrisponde un punto omologo sul piano del fotogramma (spazio immagine). Un ruolo chiave nella formazione di immagini sui sensori passivi è svolto dai sistemi ottici. Per questo motivo, è necessario fare una breve introduzione sul concetto di imaging attraverso un sistema ottico prima di illustrare i concetti di base per la costruzione di un piano di volo nel caso in cui venga utilizzata una camera di tipo “frame”. Il principio di proiezione dell'immagine, che rappresenta il processo di imaging ideale di un oggetto reale sul piano del sensore (e quindi nell'immagine), si basa esclusivamente su un principio geometrico. Lenti sferiche il cui spessore sia tanto piccolo da poter essere giudicato trascurabile rispetto alle altre grandezze in gioco (raggi di curvatura, distanze focali, ecc.) sono dette lenti sottili e sono quelle maggiormente usate nel campo fotogrammetrico. Si definisce centro ottico O di una lente sottile il punto individuato dall’intersezione della lente con l’asse ottico. Esso non provoca nessuna deviazione che non provochi alcuna deviazione a qualunque raggio luminoso che lo intercetti. Una lente sottile può essere rappresentata convenzionalmente con un segmento, perpendicolare all’asse ottico, il cui punto di mezzo centrale rappresenta il centro ottico. Quando una lente convergente intercetta un fascio di raggi luminosi, con direzione parallela al suo asse ottico e provenienti dalla parte sinistra della lente stessa, questi emergono dalla parte opposta della lente formando un cono luminoso convergente in un punto sull’asse ottico della lente chiamato fuoco (secondo fuoco). Dal momento che il fascio di raggi paralleli all’asse ottico può giungere sulla lente dalle due parti opposte, ne consegue che ogni lente possiede due fuochi, indicati con F1 e F2 (primo e secondo fuoco), entrambi sull’asse ottico, ma dalle parti opposte della lente. In una lente sottile, anche con CORSO DI GIS E 3D CITY MODEL 61  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città raggi di curvatura diversi, la distanza tra i due fuochi e il centro della lente è uguale, è chiamata distanza focale e viene indicata con f. Considerando una sorgente luminosa puntiforme P (Figura 4.1.), alla sinistra del primo fuoco di una lente convergente, la sua immagine (reale) sarà il punto P’ dove si intersecano i raggi provenienti da P dopo essere stati rifratti attraverso la lente. Per definire la posizione dell’immagine P’ si possono considerare almeno due dei tre seguenti raggi luminosi, scelti tra gli infiniti che escono da P: - il raggio AB parallelo all’asse ottico (1), che è deviato dalla lente in modo da passare per il secondo fuoco F2 della lente; - il raggio AF1 diretto sul primo fuoco della lente (2), che, intercettando la lente, viene deviato per emergere parallelo all’asse ottico; - il raggio AO che, attraversando il centro ottico della lente (3) prosegue senza subire alcuna deviazione. La costruzione delle immagini nelle lenti sottili viene, quindi, confermate anche dall’esperienza. Semplificando il presupposto di Gauss e prendendo in considerazione lenti sottili (come quelle utilizzate nelle fotocamere per la fotogrammetria), l'equazione di Huygens diventa (Kraus 2007): 1 1 1   [4.1] c d0 di dove c è la lunghezza focale, d0 è la distanza tra l'oggetto e il centro dell'obiettivo e di è la distanza tra l'immagine dell'oggetto e il centro dell'obiettivo (Figura 4.1). Figura 4.1 – Processo di formazione dell’immagine nelle lenti sottili CORSO DI GIS E 3D CITY MODEL 62  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città Nel caso d0 sia molto grande rispetto a c e di, il piano in cui l'immagine è formata è praticamente coincidente con il piano focale dell'obiettivo (c≈di). Lo schema proiettivo che si realizza è quindi una prospettiva centrale. Pertanto, nel caso di terreno perfettamente piatto e orizzontale, l'asse della fotocamera sarà in posizione nadir (o verticale) e la scala dell'immagine sarà dedotta dalla relazione: Z mb  [4.2] c essendo Z l'altezza di volo dal suolo. In generale, gli obiettivi fotografici reali non rispecchiano la condizione geometrica ideale. Infatti, come mostrato in figura 1.8, il raggio che proviene dall’oggetto incontra il primo punto nodale formando un angolo α esterno all’obiettivo e diverso dall’angolo interno α′ che il raggio uscente dal secondo punto nodale forma con l’asse ottico. Il punto immagine I’ di I si forma in I’’ a causa della distorsione introdotta dall’obiettivo, distorsione che non è solo nel piano della figura 4.2. Figura 4.2 – Processo di formazione reale dell’immagine nelle lenti sottili La distorsione δI può essere sempre scomposta in due componenti: una lungo la congiungente P con I’, detta radiale δr, e una normale a questa direzione, detta tangenziale δt. La distorsione radiale è prevalente rispetto a quella tangenziale (nei buoni obiettivi si può ritenere quest’ultima nell’ordine del 5% rispetto all’altra). La distorsione radiale può essere di due tipi: a barilotto o a cuscinetto, come mostrato in figura 1.10. CORSO DI GIS E 3D CITY MODEL 63  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città Figura 4.3 – Distorsione radiale: immagine corretta (sinistra), distorsione a barilotto (centro) e a cuscinetto (destra) 4.2.2. IL FOTOGRAMMA Il termine fotogramma si riferisce all’immagine realizzata con una camera metrica. In una camera fotogrammetrica analogica, l’immagine osservata viene proiettata sul piano focale dove è situata la pellicola; la luce riflessa della scena ripresa sensibilizza la pellicola, ovvero le cariche elettriche che costituiscono la luce, (i fotoni) colpiscono i granuli di alogenuro d’argento dell’emulsione fotografica e li modificano in modo tale che quando la pellicola subisce il procedimento chimico di sviluppo, essi si anneriscono tanto più quanto è il numero di fotoni che li ha colpiti. Nel momento della presa si forma nell’emulsione della pellicola l’immagine latente che, costituisce, una discretizzazione dell’immagine della realtà secondo il tessuto irregolare dei granuli di alogenuro d’argento. Lo sviluppo delle immagini è quello tipico delle pellicole. Nei fotogrammi delle fotocamere a pellicola era solito riportare una strip in cui si riportavano i valori della focale, nonché l’altezza di volo, la data, l’ora, la posizione della bolla sferica, i numeri di serie e la progressione dei fotogrammi scattati. Su ciascun fotogramma sono presenti 4 reperes al centro di ogni lato e 4 marche fiduciali in corrispondenza degli spigoli; collegando le due coppie di reperes contrapposti e/o le due coppie di marche fiduciali contrapposte, le linee si intersecano nel centro del fotogramma, ovvero nel punto principale. Con l’avvento dei fotogrammi digitali sono scomparsi tutti i riferimenti che erano presenti sula strip sono scomparsi e, tutte le informazioni, sono contenute nell’EXIF file associato a ciascun fotogramma digitale. Inoltre, prima dell’avvento delle camere digitali, il formato dei fotogrammi delle camere per presa aerea era standardizzato: dimensioni del fotogramma 230 x 230 mm e lunghezza focale di 150 mm (fotocamera grandangolare). Talvolta veniva usata una camera aerofotogrammetrica ad angolo normale (300 mm) quando era necessario volare ad un'altitudine relativa più elevata o nel caso della produzione ortofoto per mitigare le distorsioni prospettiche. L’impatto della scelta di un obiettivo condiziona l’angolo di vista (denominato Field of View –FOV): più la lunghezza focale CORSO DI GIS E 3D CITY MODEL 64  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città è corta più è ampia la scena, più la focale è lunga meno area della scena osserviamo, come mostrato di seguito (figura 4.4). Figura 4.4 –Variazione dell’angolo di vista in relazione alla focale Mentre nelle vecchie fotocamere analogiche la conoscenza dell'obiettivo focale era sufficiente per fare un piano di volo, nelle moderne fotocamere aeree digitali (DAC) la combinazione di lunghezza focale, dimensione dei pixel (di solito quadrata) e il formato del sensore determina il suo profilo operativo. In effetti, nel DAC questi parametri possono differire da un prototipo di telecamera a un altro. Il parametro più importante è ora la distanza di campionamento del terreno (GSD), che viene calcolata con la formula (Pepe et al., 2018): Z GSD   CCD pixel size [4.3] c 4.2.3. EQUAZIONE DI COLLINEARITÀ 0 𝑥0 𝑦0 Sia ℝ2 ∖ {0} il piano privato dell’origine 𝑂 = ( ). Due vettori 𝑋 = (𝑥 ) e 𝑌 = (𝑦 ) in ℝ2 ∖ 0 1 1 {0} si dicono collineari se appartengono alla stessa retta per l’origine, ossia se sono uno multiplo dell’altro per uno scalare non nullo. L’insieme dei raggi luminosi che, passando attraverso l’obiettivo, forma l’immagine fotografica costituisce un fascio la cui forma dipende dalle caratteristiche geometriche della fotocamera utilizzata. Tenendo conto dello schema indicato in figura 4.5, per poter effettuare delle misure è indispensabile conoscere i parametri di orientamento interno. CORSO DI GIS E 3D CITY MODEL 65  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città Figura 4.5 –Schema di presa ed orientamento interno Consideriamo un fotogramma inquadrato in un sistema di riferimento 𝑋, 𝑌, 𝑍 mente con 𝑥, 𝑦, 𝑧 sono indicati gli assi legati al centro di proiezione; gli assi 𝑥 e 𝑦 sono orientati rispettivamente come 𝜉 e 𝜂 mentre il punto O rappresenta l’origine. Consideriamo inoltre, un nuovo sistema di riferimento 𝑥 ∗ , 𝑦 ∗ , 𝑧 ∗ ruotato rispetto al precedente in maniera tale che gli assi di questa nuova terna risultino paralleli al sistema di riferimento 𝑋, 𝑌, 𝑍, ovvero: 𝑥∗ 𝑟11 𝑟12 𝑟13 𝑥 [𝑦 ] = |𝑟21 ∗ 𝑟22 𝑟23 ] ∙ [𝑦] [4.4] 𝑧∗ 𝑟31 𝑟32 𝑟33 𝑧 Un generico punto P nello spazio di coordinate 𝑋𝑃 , 𝑌𝑃 , 𝑍𝑃 , sul fotogramma avrà coordinate 𝑥𝑃 , 𝑦𝑃 , 𝑧𝑃. Nell’ipotesi che il punto principale PP coincide il centro fiduciale FC; di conseguenza, le coordinate 𝜉0 e 𝜂0 sono nulle. CORSO DI GIS E 3D CITY MODEL 66  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città Figura 4.6 – Rappresentazione dello schema di presa per la determinazione delle equazioni di collinearità Dalla similitudine dei triangoli, è possibile ricavare le seguenti espressioni: 𝑋 − 𝑋0 𝑥 ∗ = [4.5] 𝑍 − 𝑍0 𝑧 ∗ 𝑌 − 𝑌0 𝑦 ∗ = [4.6] 𝑍 − 𝑍0 𝑧 ∗ Esplicitando la 2 e la 3 si ottiene: 𝑥∗ X = 𝑋0 + (𝑍 − 𝑍0 ) ∗ { 𝑧 [4.7] 𝑥∗ 𝑌 = 𝑌0 + (𝑍 − 𝑍0 ) ∗ 𝑧 Nel sistema di riferimento 𝑥, 𝑦, 𝑧 i punti immagine positiva avranno le seguenti coordinate: CORSO DI GIS E 3D CITY MODEL 67  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città 𝑥=𝜉 {𝑦=𝜂 [4.8] 𝑧 = −𝑐 Sostituendo la 5 nella 4 e tenendo conto della relazione 1, è possibile ottenere una nuova relazione in funzione di 𝑋 e 𝑌: 𝑟11 𝜉 + 𝑟12 𝜂 − 𝑟13 𝑐 X = 𝑋0 + (𝑍 − 𝑍0 ) 𝑟31 𝜉 + 𝑟32 𝜂 − 𝑟33 𝑐 [4.9] 𝑟21 𝜉 + 𝑟22 𝜂 − 𝑟32 𝑐 Y = 𝑌0 + (𝑍 − 𝑍0 ) { 𝑟31 𝜉 + 𝑟32 𝜂 − 𝑟33 𝑐 Generalizzando al caso in cui 𝜉0 e 𝜂0 sono non nulle, si ha: 𝑟11 (𝜉 − 𝜉0 ) + 𝑟12 (𝜂 − 𝜂0 ) − 𝑟13 𝑐 𝑋 = 𝑋0 + (𝑍 − 𝑍0 ) 𝑟31 (𝜉 − 𝜉0 ) + 𝑟32 (𝜂 − 𝜂0 ) − 𝑟33 𝑐 [4.10] 𝑟21 (𝜉 − 𝜉0 ) + 𝑟22 (𝜂 − 𝜂0 ) − 𝑟32 𝑐 Y = 𝑌0 + (𝑍 − 𝑍0 ) { 𝑟31 (𝜉 − 𝜉0 ) + 𝑟32 (𝜂 − 𝜂0 ) − 𝑟33 𝑐 Infine, riportando le incognite al caso 𝜉 e 𝜂, si ha: 𝑟11 (𝑋 − 𝑋0 ) + 𝑟21 (𝑌 − 𝑌0 ) + 𝑟31 (𝑍 − 𝑍0 ) 𝜉 = 𝜉0 − 𝑐 𝑟13 (𝑋 − 𝑋0 ) + 𝑟23 (𝑌 − 𝑌0 ) + 𝑟33 (𝑍 − 𝑍0 ) [4.11] 𝑟12 (𝑋 − 𝑋0 ) + 𝑟22 (𝑌 − 𝑌0 ) + 𝑟32 (𝑍 − 𝑍0 ) 𝜂 = 𝜂0 − 𝑐 { 𝑟31 (𝑋 − 𝑋0 ) + 𝑟23 (𝑌 − 𝑌0 ) + 𝑟33 (𝑍 − 𝑍0 ) Quest’ultima equazione mostra che a causa della presenza della coordinata 𝑍 al secondo membro, per ogni punto immagine esistono infiniti possibili punti oggetto. Di conseguenza, è impossibile ricostruire la geometria spaziale di un oggetto a partire da un solo fotogramma. È necessario quindi di disporre di un ulteriore fotogramma che osserva lo stesso oggetto. 4.2.4. PIANIFICAZIONE DI VOLO TRADIZIONALE PER LA STEREO- FOTOGRAMMETRIA La stereoscopia è la tecnica per la riproduzione degli effetti tridimensionali a partire da almeno due fotogrammi (piani) che presentano una predeterminata sovrapposizione. La costruzione della nuvola di punti si sviluppa per fasi operative distinte e successive: - Copertura aerofotogrammetrica; CORSO DI GIS E 3D CITY MODEL 68  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città - Inquadramento geometrico generale; - Restituzione fotogrammetrica numerica; Lo schema della ripresa aerea è mostrato in figura 4.7: in esso l’aereo segue una traiettoria rettilinea a velocità costante, a una certa altezza media H dal terreno; all’istante t1 il centro di presa della camera si trova in O1 e riprende il primo fotogramma, e all’istante t2 (dopo aver percorso un certo spazio che costituisce la base di presa B) esso si troverà nella posizione O2 dove riprenderà il secondo fotogramma (Cannarozzo et al., 2012). Figura 4.7 – Esempio di schema di presa L’aereo sorvola l’area di interesse ripetendo lo schema precedente ed eseguendo le prese nell’ambito di una sequenza di percorsi rettilinei affiancati lungo direzioni parallele secondo lo schema illustrato in figura 4.8. In effetti, più schemi si possono realizzare per coprire di fotogrammi l’area da cartografare. CORSO DI GIS E 3D CITY MODEL 69  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città Figura 4.8 –Possibili schemi di acquisizione per missioni di fotogrammetria aerea: “zig-zag” (a); “jump line” (b); “one way” or “racetrack” (c). I fotogrammi ripresi lungo uno stesso percorso rettilineo costituiscono una strisciata. L’insieme di più strisciate viene chiamato blocco fotogrammetrico. Oltre al ricoprimento longitudinale (overlap μ (in generale il 60÷65%), è necessario che l’interasse i, tra due strisciate adiacenti, venga stabilito in modo che ci sia anche un ricoprimento laterale ε (sidelap) con valore compreso nell’intervallo del 15%÷30%) dell’abbracciamento L del fotogramma (nell’ipotesi di fotogramma di forma quadrata), per evitare di avere «buchi» nella copertura del territorio (figure 4.9). Le prese fotografiche del terreno devono avvenire in modo tale che tutta l’area da rilevare rimanga scomposta in modelli stereoscopici tali che ogni punto del terreno compaia, come detto, almeno su due fotogrammi. Figura 4.9 –Ricoprimento trasversale e longitudinale CORSO DI GIS E 3D CITY MODEL 70  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città 4.2.5. LA GESTIONE DEL VOLO FOTOGRAMMETRICO L’obiettivo primario del piano di volo è identificare i centri di presa di ogni fotogramma al fine di rispettare i valori fotogrammetrici di overlap e sidelap. Per fare questo, la camera fotogrammetrica è legata ad un sistema di navigazione GNSS e ad un software di gestione del piano di volo. Un osservatore, posto a bordo dell’aereo accanto al pilota, verificherà il rispetto del piano di volo. Figura 4.10 –Controller montato a bordo dell’aereo per gestione del piano di volo Una problematica che potrebbe insorgere nella missione fotogrammetrica è dovuto alla presenza di piccole nubi, incendi o foschie che potrebbero inficiare sulla qualità dei fotogrammi. Altre problematiche che potrebbero insorgere riguardano i movimenti impressi all’aereo dai venti. Infatti, per effetto del vento laterale che preme sull’impennaggio verticale della fusoliera, l’aereo non si muove secondo la direzione dell’asse della fusoliera, ma nella direzione della risultante R tra l’azione del vento Vv e la spinta di propulsione dell’aereo Vp. Nel nostro contesto la deriva è l’angolo δ, compreso tra la direzione della risultante R (cioè la rotta dell’aereo) e l’asse della fusoliera dell’aereo. Se la camera fosse tenuta in posizione fissa rispetto alla fusoliera, si otterrebbe una strisciata con fotogrammi male orientati, come quella rappresentata in figura 4.11a, in cui il ricoprimento longitudinale sarebbe errato. Grazie ad una piattaforma inerziale che consente la CORSO DI GIS E 3D CITY MODEL 71  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città rotazione all’asse verticale della camera, i fotogrammi possono essere acquisiti nel modo corretto (figura 4.11b). (a) (b) Figura 4.11 –Effetto della deriva sui fotogrammi della strisciata: senza la correzione della camera (a); con correzione ella camera (b). 4.2.6. INQUADRAMENTO GEOMETRICO GENERALE Per svolgere le operazioni di orientamento esterno delle immagini digitali primarie è necessario determinare un insieme di punti d’appoggio, denominati GCPs (Ground Control Points) e di controllo, o anche chiamati CPs (Check Points). I risultati dell’orientamento esterno dipendono anche dalla precisione con cui tali punti sono stati determinati. In generale i GCP e i CP devono essere determinati mediante adeguate operazioni di rilievo topografico e geodetico sul terreno. Tali tecniche sono ormai essenzialmente svolte con tecniche GNSS (Global Navigation Satellite System) in varie modalità (rilievo rapido-statico, RTK, ecc.) spesso appoggiate a reti di stazioni permanenti. CORSO DI GIS E 3D CITY MODEL 72  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città Figura 4.12 – Segnalazione e procedura di rilevamento dei GCPs 4.2.7. TRIANGOLAZIONE AEREA Per triangolazione aerea (TA) si intende quel procedimento analitico che consente l’orientamento esterno simultaneo di un blocco di fotogrammi, vincolando le posizioni relative grazie a punti di legame. I punti di legame (detti anche Tie points) sono punti comuni a due o più fotogrammi e di coordinate oggetto incognite. L’orientamento esterno con la triangolazione aerea consente di ridurre drasticamente il numero di punti fotogrammetrici di appoggio necessari. Il metodo di triangolazione aerea più diffuso è il Bundle Adjustment, ovvero metodo a stelle proiettive. Figura 4.13. Principio della compensazione di un blocco a stelle proiettive Le stelle proiettive vengono traslate e ruotate, in modo che i raggi si intersechino al meglio in corrispondenza dei punti di legame e passino il più possibile per i punti di appoggio (Kraus, 1994). Il sistema di equazioni alle osservazioni presenta, come incognite, i 6 parametri di orientamento esterno per ogni fotogramma e le tre coordinate oggetto di ogni Tie point. Il sistema, ridondante, si risolve con il metodo della compensazione ai minimi quadrati che fornisce le correzioni da CORSO DI GIS E 3D CITY MODEL 73  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città apportare ai valori approssimati delle incognite e i valori dei residui sulle osservazioni. A tale scopo le equazioni devono essere linearizzate. 4.3. Le fotocamere oblique (frame cameras) Negli ultimi anni, nel settore aerofotogrammetrico, si sono sviluppate tecnologie costituite da multi camere digitali oblique, in grado non solo di riprendere lo stesso target (edifici, antenne e pali) da diverse angolazioni, ma di determinarne, grazie ad opportuni software dedicati, anche la geometria. In generale, il sistema aviotrasportato è costituito (Figura 4.14) da quattro camere digitali metriche di tipo “frame” montate in posizione obliqua e disposte secondo un orientamento Nord, Est, Sud, Ovest nonché da un’ulteriore camera in posizione nadirale. L’integrazione delle camere aerofotogrammetriche con sistemi GNSS-INS consente la georeferenziazione diretta delle immagini acquisite. La georeferenziazione diretta permette di conoscere i 6 parametri di orientamento esterno (posizione ed assetto) del singolo fotogramma in uno specifico sistema di riferimento. Figura 4.14. Osservazione di un target da più viste CORSO DI GIS E 3D CITY MODEL 74  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città Un esempio di sistema di fotocamera obliqua è il Leica RCD30 Oblique Penta System, che consiste di cinque fotocamere: quattro fotocamere inclinate a 35° e lungo quattro direzioni ortogonali, e una fotocamera looking-nadir. Altri esempi di sistema di fotocamera obliqua sono MIDAS della società TrackAir (Madani, 2012), IGI DigiCAM Penta (Jacobsen & Gerke, 2012), e Pictometry Vexcel Osprey (Gruber & Walker, 2014). La forma della copertura del terreno catturata simultaneamente dalle cinque fotocamere (figura 4.15) sembra una "Croce di Malta" (Petrie, 2009). Figura 4.15. Leica RCD30 Oblique Penta System In tale configurazione, la formula per il calcolo della GSD è riportata in Höhle (2008; 2111) e qui richiamata: Z cos(  t ) GSD  pixelsize'  [4.12] c cos  dove  è l'angolo tra la linea diretta dall'obiettivo a un bersaglio e la verticale, mentre t è l'angolo di inclinazione della telecamera considerata. Figura 4.16. Schema di acquisizione con fotocamere oblique CORSO DI GIS E 3D CITY MODEL 75  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città 4.3.1. Impiego dell’approccio Structure-from-Motion (SfM) alle frame-camera L'approccio Structure from Motion (SfM) è diventato molto popolare nel campo fotogrammetrcio grazie alla capacità di determinare i parametri di orientamento esterno senza alcuna conoscenza a priori delle posizioni approssimative delle fotocamere e dei punti 3D. La tecnica SfM richiede, per realizzare modelli 3D, un blocco di immagini con un alto grado di sovrapposizione che catturi l'intera struttura 3D della scena vista da più posizioni. Figura 4.17. Principio di acquisizione in ambiente SfM Al fine di aumentare la densità della nuvola di punti generata dall’approccio SfM, diversi algoritmi denominati Multi View Stereo (MVS) sono stati sviluppati (Skarlatos e Kiparissi 2012). A causa delle grandi dimensioni delle immagini aeree digitali generate dai sensori aerei e, di conseguenza, dei problemi di computazione, queste immagini sono state raramente utilizzate per la costruzione di modelli 3D attraverso l'approccio SfM-MVS. Inoltre, i valori di sovrapposizione utilizzati nei tradizionali rilievi aerei porterebbero ad una geometria di acquisizione piuttosto debole. Per questi motivi, l'uso di immagini aeree generate da sensori aerei nell'approccio SfM- MVS è stato limitato. Tuttavia, tenendo conto della possibilità di acquisire nadir aerei e immagini oblique secondo la multi-vista, dell'aumento del calcolo ad alte prestazioni e dell'uso della CORSO DI GIS E 3D CITY MODEL 76  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città georeferenziazione diretta, Pepe et al. 2016 hanno costruito un modello 3D della città svizzera di Nöllen (figura 4.18). Figura 4.18. Posizioni delle fotocamere su modello di elevazione. Pertanto, i diversi passaggi necessari per ottenere il modello 3D possono essere riassunti come segue: - Allineamento delle immagini (Alignment of the images). Questo compito permette di determinare la posizione e l'orientamento di ogni immagine. Inoltre, il software SfM produce una rada nuvola di punti. I due passaggi principali che portano alla ricostruzione del modello possono essere riassunti come segue (Snavely 2008): o Ricerca della corrispondenza:  Estrazione delle caratteristiche (Feature extraction): questo compito può essere realizzato, ad esempio, utilizzando il rilevatore di caratteristiche di trasformazione delle caratteristiche in scala (SIFT) (Lowe 2004) che è un algoritmo di rilevamento delle caratteristiche nella visione artificiale per rilevare e descrivere le caratteristiche locali nelle immagini.  Corrispondenza (Matching): il sistema cerca le corrispondenze delle caratteristiche trovando la caratteristica più simile nell'immagine Ia per ogni caratteristica nell'immagine Ib, utilizzando una metrica di somiglianza che ne confronta l'aspetto. CORSO DI GIS E 3D CITY MODEL 77  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città  Verifica geometrica (Geometric verification): questa fase verifica le coppie di immagini potenzialmente sovrapposte utilizzando una tecnica di stima robusta, come il RANSAC (RANdomSAmple Consensus) (Fischler e Bolles 1987). o Fase di ricostruzione;  Inizializzazione (Initialization): il sistema inizia con una ricostruzione a due viste accuratamente selezionata.  Registrazione delle immagini (Image registration): in questo compito, utilizzando le corrispondenze delle caratteristiche, si possono registrare nuove immagini al modello corrente risolvendo il problema del Perspective-n-Point (PnP) (Fischler e Bolles 1987).  Triangolazione (Triangulation): una nuova immagine registrata deve osservare i punti di scena esistenti e, anche in questa fase, l'uso di una tecnica di stima robusta è fondamentale.  Bundle Adjustment (BA): è un metodo di perfezionamento utilizzato per migliorare le soluzioni SfM  Costruire nuvole dense (dense point cloud). Sulla base delle posizioni stimate della telecamera, è possibile calcolare le informazioni di profondità per ogni immagine in modo da poterle combinare in un'unica nuvola di punti densi. Sulla base dell'algoritmo MVS, questo compito permette di aumentare la densità della nuvola di punti generata nel processo SfM.  Costruire la mesh. Il processo di trasformazione dalla nuvola di punti in una mesh triangolare a superficie singola può essere realizzato con diversi approcci algoritmici in relazione alla morfologia dell'oggetto o dell'area in esame. 4.4. La georeferenziazione diretta (Direct Georeferencing – DG) 4.4.1. Principi L’integrazione dei sistemi di posizionamento satellitare (GNSS) con quelli inerziali (INS) è di rilevante importanza in diversi campi della geomatica come, ad esempio, in aerofotogrammetria e nel Remote Sensing, in quanto la tecnica GNSS/INS consente di ricostruire, istante per istante, la posizione e l’assetto del sensore (figura 4.19) CORSO DI GIS E 3D CITY MODEL 78  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città Figura 4.19- I due differenti sistemi angolari: (R, P, H) e (, , k) Il punto di forza dell’integrazione GNSS/INS è il loro comportamento complementare: il GNSS aiuta il sistema di navigazione inerziale a fornire con elevata accuratezza la posizione mentre, dall’altro lato, quello inerziale integra il sistema GNSS fornendo un’accurata posizione iniziale ed informazione di velocità quando si verifica una perdita del segnale proveniente dai satelliti, continuando ad avere informazioni sulla traiettoria. PROPRIETÀ INS GPS Frequenza output Alta (50-100 Hz) Bassa (1-10 Hz) Accuratezza a lungo termine Bassa Alta Accuratezza a lungo termine Alta Alta Dipendente dalla diretta Indipendenza Unità autonoma osservabilità dei satelliti Capacità di misurare Preciso e affidabile Preciso e affidabile l’assetto Tabella 3.1 – Confronto vantaggi e svantaggi tra INS e GPS CORSO DI GIS E 3D CITY MODEL 79  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città 4.4.2. La navigazione inerziale La navigazione inerziale si basa sulle leggi della dinamica di Newton: dal secondo principio, considerato un corpo di massa m, sottoposto ad una forza esterna 𝐹̅ , il suo moto sarà definito, in ogni istante, da un vettore accelerazione 𝑎̅ parallelo ad 𝐹̅ , dato dalla relazione: 𝑎̅ = 𝐹̅ ⁄𝑚 [4.13] Nell’ipotesi di un modello di Terra sferica, con asse di rotazione inerziale (trascurando il moto di rivoluzione della Terra intorno al Sole e il moto di precessione dell’asse terrestre) e con riferimento ad una terna, i cui gli assi x-y giacciono sul piano equatoriale, fissa rispetto alle stelle e l’asse z coincide con l’asse terrestre, la formula [4.13] diviene (Nastro, 2004): 𝑑 2 𝑅̅ [4.14] 𝑎̅ = ( ) = 𝑓 ̅ + 𝐺̅ 𝑑𝑡 2 𝑖 che rappresenta l’equazione fondamentale della navigazione inerziale, in cui: 𝑅̅ vettore posizione del mobile in un generico punto; 𝑎̅ vettore accelerazione assoluta o inerziale, riferita alla terna inerziale a cui è sottoposta una massa m posta in un generico punto; 𝑓̅ vettore somma di tutte le accelerazioni applicate alla massa m; 𝐺̅ vettore accelerazione gravitazionale. Per ricavare il vettore posizione 𝑅̅ occorrerebbe misurare l’accelerazione assoluta 𝑎̅; integrando detta misura, tenendo conto del valore della velocità iniziale, si ricaverebbe la velocità; integrando ancora la velocità calcolata, nota la posizione iniziale, si ricaverebbe la posizione corrente. Purtroppo, per il principio di equivalenza di Einstein, è impossibile distinguere tra gli effetti dell’accelerazione inerziale e quella del campo gravitazionale e, di conseguenza, un qualsiasi misuratore di accelerazioni sarà influenzato dal campo gravitazionale 𝐺̅. Pertanto è necessario sommare alla misura effettuata 𝑓 ̅ il valore del campo gravitazionale 𝐺̅ , che è funzione della località. Il modulo e la direzione del vettore 𝑓 ̅ possono essere ricavati con l’ausilio di accelerometri e giroscopi. Gli accelerometri e i giroscopi sono impiegati per misurare, rispettivamente, l’accelerazione lineare e la velocità angolare (Woodman, 2007). Il principio base degli accelerometri consiste nel misurare la forza agente su una massa di prova (proof mass). Si distinguono principalmente due tipi di sensori: open-loop o closed-loop. Nel primo caso gli accelerometri misurano direttamente lo spostamento della massa di prova dovuto alle forze CORSO DI GIS E 3D CITY MODEL 80  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città esterne agenti sul sensore mentre i sensori closed-loop mantengono la massa di prova in uno stato di equilibrio generando una forza opposta a quella applicata. Questo è ottenuto attraverso un sistema di controllo retroazionato (feedback system) di tipo elettrico o magnetico. Il principale vantaggio dell’architettura closed-loop è la miglior linearità della misura rispetto alle forze applicate (Biagi, 2009). Il giroscopio è invece, uno strumento in grado di misurare la velocità angolare della rotazione impressa al sensore rispetto a un sistema di riferimento inerziale; in relazione alle sue caratteristiche, esso può essere di tipo meccanico, ottico e MEMS (Micro-Electrical Mechanical System). Un giroscopio meccanico convenzionale (Figura 4.20) è costituito da un rotore montato su due giunti cardanici e, questa configurazione, segue la rotazione intorno ai tre assi. Per la conservazione del momento angolare, la rotazione del rotore rimarrà con un orientamento globale costante poiché cambiano gli angoli tra le sospensioni cardaniche; questi angoli sono, quindi, misurati al fine di monitorare l’orientamento. Figura 4.20 – Giroscopio meccanico I giroscopi ottici, invece, utilizzando l'interferenza della luce e sfruttando l’effetto di Sagnac, misurano la velocità angolare (Woodman, 2007). In particolare, quelli attualmente utilizzati sono a fibra ottica FOG (Fiber Optic Gyroscope) i quali garantiscono un’elevata precisione (Biagi, 2009). Infine, per misurare la velocità angolare i giroscopi MEMS sfruttano l’effetto di Coriolis, misurato usando elementi geometrici vibranti. Sebbene questi sensori non possano raggiungere la stessa accuratezza di quelli ottici, hanno vari vantaggi: dimensioni e peso ridotto, basso consumo energetico e sono economici da produrre. Pertanto, nelle piattaforme inerziali, a seconda dell’accuratezza attesa, viene scelto quello ottico o MEMS. CORSO DI GIS E 3D CITY MODEL 81  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città Un sistema di navigazione inerziale strapdown utilizza accelerometri ortogonali e terne di giroscopi rigidamente fissati agli assi del veicolo in movimento. Il movimento angolare del sistema è misurato continuamente con i sensori di velocità; gli accelerometri non rimangono stabili nello spazio ma seguono il moto del veicolo. L'algoritmo di una piattaforma strapdown è rappresentato in figura 4.21. Figura 4.21 - Algoritmo di una piattaforma inerziale strapdown (Woodman, 2007) 4.4.3. Sistemi di coordinate In letteratura, sono presenti diversi sistemi di coordinate utilizzati per rappresentare gli stati dei sistemi di navigazione e i relativi parametri cinematici, di cui si riporta di seguito, una classificazione. 1) Il body frame (b-frame): il b-frame (Figura 4.22), noto anche come frame dell’aeromobile, è un sistema di assi ortogonali che ha l’origine, P, nel centro del veicolo (coincidente con quella del frame locale di navigazione) ed è allineato con l'asse di beccheggio, xb,, l’asse di rollio, yb, e l’asse di imbardata, zb, del veicolo su cui è installato il sistema di navigazione. Il body frame è essenziale nella navigazione in quanto descrive il moto dell’aeromobile; inoltre, va sottolineato che tutti i sensori strapdown misurano il moto di questo frame rispetto ad un sistema inerziale generico. CORSO DI GIS E 3D CITY MODEL 82  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città Figura 4.22 – Body frame 2) Il frame inerziale (i-frame) o ECI (Earth Centered Inertial): l'i-frame (Figura 3.10) ha l’origine nel centro della Terra e gli assi sono non rotanti rispetto allo spazio inerziale. Per la navigazione aerea, il frame ECI è un’approssimazione di quell’inerziale newtoniano. L’asse xi giace nel piano equatoriale ed è diretto verso l’equinozio di primavera, l’asse zi è allineato con quello di rotazione terrestre, infine l’asse yi completa la terna levogira. 3) Il frame terrestre (e-frame) o ECEF (Earth Centered Earth Fixed): l’e-frame (Figura 3.10) ha origine nel centro di massa della Terra e ruota con essa; l’asse ze, coincidente con l’asse di rotazione terrestre, punta verso la direzione del Polo Convenzionale Terrestre (CTP - Conventional Terrestrial Pole); l’asse xe giace nel piano equatoriale e punta verso l’intersezione dell’equatore del CTP con il meridiano di riferimento (meridiano zero definito dal BIH); infine, l’asse ye completa il sistema levogiro. Rispetto al frame inerziale l’e-frame ruota, a una velocità ωie, attorno all'asse zi. Il WGS84 è un frame ECEF di uso comune che definisce il modello della Terra e il relativo ellissoide di riferimento. In tale frame le coordinate di un punto possono essere espresse in coordinate geodetiche o cartesiane, legate dalla seguente trasformazione: 𝑥𝑒 = (𝑁𝑃𝑉 + ℎ) 𝑐𝑜𝑠𝜑 𝑐𝑜𝑠𝜆 { 𝑦𝑒 = (𝑁𝑃𝑉 + ℎ) 𝑐𝑜𝑠𝜑 𝑠𝑖𝑛𝜆 [4.15] 𝑧𝑒 = [(1 − 𝑒 2 ) 𝑁𝑃𝑉 + ℎ]𝑠𝑖𝑛𝜑 avendo indicato con 𝑁𝑃𝑉 il raggio di curvatura del primo verticale. 4) Il frame di navigazione (n-frame): l'n-frame è un frame locale, geodetico o geografico, che risulta fisso rispetto al velivolo in quanto ha l’origine, P, nella posizione di navigazione CORSO DI GIS E 3D CITY MODEL 83  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città del sistema, ossia nel suo centro di massa; l’asse zn punta verso il basso, o verso l’alto, perpendicolarmente all’ellissoide di riferimento; gli assi xn e yn, invece, giacciono in un piano tangente a tale ellissoide. Se l’asse xn punta a Nord, yn a Est e zn verso il basso, il frame prende il nome di NED (North-East-Down); se invece, zn punta verso l’alto, il frame è definito ENU (East-North- Up), come mostrato in figura 3.10 in cui 𝑥𝑛 ≡ 𝑥𝐸 ; 𝑦𝑛 ≡ 𝑦𝑁 ; 𝑧𝑛 ≡ 𝑧𝑈. 5) Il wander azimuth navigation frame (p-frame): il p-frame (Figura 4.23) viene adottato per evitare le singolarità nel calcolo che si verificano ai poli per l’n-frame. L’asse zp coincide ancora con la verticale, come per l’n-frame, gli assi xp e yp giacciono sempre nel piano orizzontale ma sono ruotati, intorno alla verticale, di un angolo 𝛼, detto angolo di deriva; se esso è nullo la terna wander coincide con quella di navigazione. La terna wander può ruotare rispetto a quella di navigazione, intorno all’asse z, con una velocità angolare 𝛼̇ , di conseguenza, l’angolo di deriva è variabile nel tempo. Questo frame non viene chiamato w- frame (w per wander) per evitare confusione tra w e ω che rappresenta la velocità di rotazione. La lettera p in p-frame è sinonimo di piattaforma, infatti il wander azimut navigation frame è di livello locale e quindi costituisce una piattaforma orizzontale (Zhang, 2012). Figura 4.23 – Rappresentazione dei diversi frame CORSO DI GIS E 3D CITY MODEL 84  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città 4.4.4. EQUAZIONE DELLA NAVIGAZIONE INERZIALE Un sistema di navigazione inerziale consente di ricavare la posizione, la velocità e l’assetto di un veicolo, rispetto a un sistema di coordinate di riferimento, attraverso la soluzione dell’equazione fondamentale [4.13]; il modo con cui si sceglie di risolvere tale equazione, è detto meccanizzazione del sistema inerziale. Se si indica con 𝑣̅ la velocità dell’aereo rispetto alla Terra, essa può definirsi come la derivata, rispetto al tempo, del vettore posizione riferito ad una terna di coordinate solidale alla Terra: 𝑑𝑅̅ 𝑣̅ = ( ) [4.16] 𝑑𝑡 𝑡 La velocità dell’aereo rispetto a una terna inerziale, per l’equazione di Coriolis, può esprimersi come: 𝑑𝑅̅ ( ) = 𝑣̅ + 𝜎̅ × 𝑅̅ [4.17] 𝑑𝑡 𝑖 dove 𝜎̅ rappresenta la velocità angolare terrestre. Nell’ipotesi che il vettore velocità angolare terrestre sia di modulo costante e resti inerzialmente fisso, derivando la precedente relazione rispetto al tempo si ha: 𝑑 2 𝑅̅ 𝑑𝑣̅ 𝑑𝑅̅ ( 2 ) = ( ) + 𝜎̅ × ( ) [4.18] 𝑑𝑡 𝑖 𝑑𝑡 𝑖 𝑑𝑡 𝑖 Gli accelerometri misurano le componenti del vettore 𝑓 ̅ rispetto ad assi vincolati alla piattaforma che ruota, rispetto allo spazio inerziale, con una velocità angolare 𝜔 ̅ somma della velocità angolare della piattaforma rispetto alla Terra, 𝜏̅ , e della velocità angolare 𝜎̅ della Terra stessa, ossia: 𝑑𝑣̅ 𝑑𝑣̅ [4.19] ( ) = ( ) + (𝜏̅ + 𝜎̅) × 𝑣̅ 𝑑𝑡 𝑖 𝑑𝑡 𝑃 Sostituendo la [4.17] e la [4.19] nell’equazione [4.18] e raggruppando i termini comuni, si ha: 𝑑 2 𝑅̅ 𝑑𝑣̅ ( 2 ) = ( ) + (𝜏̅ + 2𝜎̅) × 𝑣̅ + 𝜎̅ × (𝜎̅ × 𝑅̅) [4.20] 𝑑𝑡 𝑖 𝑑𝑡 𝑃 ove: 𝜏̅ × 𝑣̅ accelerazione centripeta dovuta al moto dell’aereo; 2𝜎̅ × 𝑣̅ accelerazione di Coriolis; 𝜎̅ × (𝜎̅ × 𝑅̅ ) accelerazione centripeta dovuta alla rotazione terrestre. CORSO DI GIS E 3D CITY MODEL 85  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città Sostituendo la [4.20] nell’equazione fondamentale della navigazione inerziale [4.14] ed esprimendo tutto rispetto ad 𝑓 ̅ si ricava: 𝑑𝑣̅ 𝑓̅ = ( ) + (𝜏̅ + 2𝜎̅) × 𝑣̅ − [𝐺̅ − 𝜎̅ × (𝜎̅ × 𝑅̅)] [4.21] 𝑑𝑡 𝑃 Riscrivendo la relazione [4.21] rispetto al punto 𝑃: 𝑑𝑣̅ ( ) = 𝑓 ̅ − (𝜏̅ + 2𝜎̅) × 𝑣̅ + 𝑔̅ [4.22] 𝑑𝑡 𝑃 Introducendo la matrice antisimmetrica 𝑨, l’equazione [4.22] può essere scritta nel seguente modo (Nastro, 2004): 𝑑𝑣̅ [4.23] ( ) = 𝑓 ̅ − 𝑨(𝜏̅ + 2𝜎̅)𝑣̅ + 𝑔̅ 𝑑𝑡 𝑃 Tale espressione, proiettata su ciascuno degli assi di un sistema di coordinate cartesiano, fornisce le seguenti relazioni: 𝑑𝑣 ( ) = 𝑓𝑥 − (𝜏𝑦 + 2𝜎𝑦 )𝑣𝑧 + (𝜏𝑧 + 2𝜎𝑧 )𝑣𝑦 + 𝑔𝑥 𝑑𝑡 𝑥 𝑑𝑣 ( ) = 𝑓𝑦 − (𝜏𝑧 + 2𝜎𝑧 )𝑣𝑥 + (𝜏𝑥 + 2𝜎𝑥 )𝑣𝑧 + 𝑔𝑦 [4.24] 𝑑𝑡 𝑦 𝑑𝑣 ( ) = 𝑓𝑧 − (𝜏𝑥 + 2𝜎𝑥 )𝑣𝑦 + (𝜏𝑦 + 2𝜎𝑦 )𝑣𝑥 + 𝑔𝑧 { 𝑑𝑡 𝑧 Le equazioni [4.24] rappresentano le componenti dell’accelerazione dell’aereo, o meglio della piattaforma inerziale, rispetto ad una terna di misura nel frame ECEF. 4.4.5. MATRICE DEI COSENI DIRETTORI Nel caso delle piattaforme strapdown cambia continuamente l’orientamento della terna di misura ad essa legata rispetto a quella esterna di riferimento, pertanto, è necessario calcolare, istante per istante, una matrice dei coseni direttori 𝑪 (DCM – Director Cosine Matrix) che permette di passare dalla terna di misura alla terna di riferimento. Un vettore 𝑉̅ , riferito ad una generica terna cartesiana (a) di assi 𝑥, 𝑦, 𝑧, può essere definito come la somma di tre vettori (Nastro, 2004): 𝑉𝑎 = 𝑉𝑥 + 𝑉𝑦 + 𝑉𝑧 = 𝑉𝑥 𝑖 + 𝑉𝑦 𝑗 + 𝑉𝑧 𝑘 [4.25] in cui ciascuno dei termini rappresenta il vettore ottenuto proiettando 𝑉𝑎 su ciascun asse della terna (a) di versori 𝑖, 𝑗, 𝑘. Lo stesso vettore, riferito ad una seconda terna (b) di assi 𝑥’, 𝑦’, 𝑧’, avente la stessa origine, può essere rappresentato da: CORSO DI GIS E 3D CITY MODEL 86  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città 𝑉𝑏 = 𝑉𝑥′ + 𝑉𝑦′ + 𝑉𝑧′ = 𝑉𝑥′ 𝑖′ + 𝑉𝑦′ 𝑗′ + 𝑉𝑧 ′ 𝑘′ [4.26] dove ciascun termine può essere ottenuto proiettando 𝑉𝑥 , 𝑉𝑦 , 𝑉𝑧 rispettivamente sugli assi 𝑥’, 𝑦’, 𝑧’ della terna (b). Pertanto è possibile scrivere: 𝑉𝑥′ = [𝑉𝑥 cos(𝑥 ′ 𝑥) ; 𝑉𝑦 cos(𝑥 ′ 𝑦) ; 𝑉𝑧 (𝑥 ′ 𝑧)] { 𝑉𝑦′ = [𝑉𝑦 cos(𝑦 ′ 𝑥) ; 𝑉𝑦 cos(𝑦 ′ 𝑦) ; 𝑉𝑧 (𝑦 ′ 𝑧)] [4.27] 𝑉𝑧′ = [𝑉𝑧 cos(𝑧 ′ 𝑥) ; 𝑉𝑦 cos(𝑧 ′ 𝑦) ; 𝑉𝑧 (𝑧 ′ 𝑧)] dove il generico cos(𝑥 ′ 𝑦) rappresenta il coseno dell’angolo che l’asse 𝑥 ′ della terna (b) forma con l’asse 𝑦 della terna (a) ed è indicato come “coseno direttore”. Le precedenti relazioni possono essere riscritte come: 𝑉𝑥′ = [𝑉𝑥 C11 ; 𝑉𝑦 C12 ; 𝑉𝑧 C13 ] { 𝑉𝑦′ = [𝑉𝑦 C21 ; 𝑉𝑦 C22 ; 𝑉𝑧 C23 ] [4.28] 𝑉𝑧′ = [𝑉𝑧 C31 ; 𝑉𝑦 C32 ; 𝑉𝑧 C33 ] che può, sinteticamente, essere espressa in forma matriciale come: 𝑉𝑏 = 𝑪𝑏𝑎 𝑉𝑎 [4.29] con 𝑪𝑏𝑎 matrice dei coseni direttori, le cui righe sono nel frame di arrivo (b), mentre le colonne sono in quello di partenza (a). La stessa matrice, che consente il passaggio da una terna (a) ad una terna (b), può essere ottenuta con tre successive rotazioni 𝜃1 , 𝜃2 , 𝜃3 , rispetto a ciascun asse. La rotazione 𝑪𝑧 (𝜃3 ) intorno all’asse z di un angolo 𝜃3 lascia invariata la direzione dell’asse z mentre gli assi x, y ruotano di un angolo 𝜃3 , ossia: 𝑐𝑜𝑠𝜃3 −𝑠𝑖𝑛𝜃3 0 𝑪𝑧 (𝜃3 ) = (−𝑠𝑖𝑛𝜃3 𝑐𝑜𝑠𝜃3 0) [4.30] 0 0 1 invece, la rotazione 𝑪𝑦 (𝜃2 ) intorno all’asse y di un angolo 𝜃2 è: 𝑐𝑜𝑠𝜃2 0 −𝑠𝑖𝑛𝜃2 𝑪𝑦 (𝜃2 ) = ( 0 1 0 ) [4.31] 𝑠𝑖𝑛𝜃2 0 𝑐𝑜𝑠𝜃2 infine, la rotazione 𝑪𝑥 (𝜃1 ) intorno all’asse x di un angolo 𝜃1 è: 1 0 0 𝑪𝑥 (𝜃1 ) = (0 𝑐𝑜𝑠𝜃1 𝑠𝑖𝑛𝜃1 ) [4.32] 0 −𝑠𝑖𝑛𝜃1 𝑐𝑜𝑠𝜃1 Nell’ipotesi che le rotazioni avvengono nell’ordine 𝜃1 , 𝜃2 , 𝜃3 , si ha: 𝑪𝑏𝑎 = 𝐶𝑧 (𝜃3 )𝐶𝑦 (𝜃2 )𝐶𝑥 (𝜃1 ) [4.33] CORSO DI GIS E 3D CITY MODEL 87  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città Sostituendo le equazioni [4.26], [4.27] e [4.28] nella [4.43] si ottiene la seguente espressione (Nastro, 2004): 𝑐𝑜𝑠𝜃2 𝑐𝑜𝑠𝜃3 𝑐𝑜𝑠𝜃1 𝑠𝑖𝑛𝜃3 + 𝑠𝑖𝑛𝜃1 𝑠𝑖𝑛𝜃2 𝑐𝑜𝑠𝜃3 𝑠𝑖𝑛𝜃1 𝑠𝑖𝑛𝜃3 − 𝑐𝑜𝑠𝜃1 𝑠𝑖𝑛𝜃2 𝑐𝑜𝑠𝜃3 𝑏 𝑪 = (− 𝑐𝑜𝑠𝜃2 𝑠𝑖𝑛𝜃3 𝑎 𝑐𝑜𝑠𝜃1 𝑐𝑜𝑠𝜃3 − 𝑠𝑖𝑛𝜃1 𝑠𝑖𝑛𝜃2 𝑠𝑖𝑛𝜃3 𝑠𝑖𝑛𝜃1 𝑐𝑜𝑠𝜃3 + 𝑐𝑜𝑠𝜃1 𝑠𝑖𝑛𝜃2 𝑠𝑖𝑛𝜃3 ) [4.44] 𝑠𝑖𝑛𝜃2 −𝑠𝑖𝑛𝜃1 𝑐𝑜𝑠𝜃2 𝑐𝑜𝑠𝜃1 𝑐𝑜𝑠𝜃2 che rappresenta la DCM per tre successive rotazioni. La matrice [4.44] è particolarmente utile in quanto, nei sistemi strapdown, consente di passare da una terna di riferimento ad un’altra. 4.4.6. ANGOLI DI EULERO In un sistema strapdown gli assi, lungo i quali sono disposti gli accelerometri coincidono con quelli dello stesso aeromobile, ovvero con il body frame (b). Gli accelerometri misurano le componenti del vettore accelerazione 𝑓𝑏̅ rispetto ad una terna solidale all’aeromobile e quindi misurano direttamente le accelerazioni alle quali è sottoposto l’aeromobile. Le uscite degli accelerometri non possono essere utilizzate direttamente ma devono essere prima trasformate in componenti lungo gli assi di una terna di riferimento (a) che dipende dalla meccanizzazione prescelta. A tal fine, le uscite dei giroscopi vengono utilizzate per ricavare, istante per istante, l’orientamento degli assi della terna (b) rispetto a quella di riferimento (a) tramite la matrice dei coseni direttori 𝑪𝑏𝑎. 𝑏 Tale matrice è funzione della velocità angolare 𝜔𝑎𝑏 (le cui componenti lungo gli assi dell’aeromobile sono 𝑝’, 𝑞’, 𝑟’ ) con la quale la terna legata al velivolo ruota rispetto alla terna di 𝑏 riferimento che è a sua volta, legata alla velocità 𝜔𝑖𝑏 (di componenti 𝑝, 𝑞, 𝑟) misurata dai giroscopi attraverso la relazione: 𝑏 𝑏 𝑏 𝜔𝑎𝑏 = 𝜔𝑖𝑏 − 𝜔𝑖𝑎 [4.45] 𝑏 dove 𝜔𝑖𝑎 (le cui componenti sono 𝑝𝑐 , 𝑞𝑐 , 𝑟𝑐 ) è la velocità angolare con la quale la terna di riferimento ruota rispetto alla terna inerziale che è pari a (Nastro, 2004): 𝑏 𝑎 𝜔𝑖𝑎 = 𝜔𝑖𝑎 𝐶𝑎𝑏 = (𝜏̅ + 𝜎̅)𝐶𝑎𝑏 [4.46] L’orientamento del body frame rispetto al frame di navigazione, che nei sistemi strapdown si identifica con la terna NED, è definito dalla conoscenza degli angoli di Eulero. Gli angoli di Eulero (Figura 3.11) si definiscono a partire dalla posizione canonica dell’aeromobile, in maniera che la terna body ruoti fino a coincidere con la terna NED (velocità di rotazioni antiorarie positive). Per passare da tale posizione a un assetto qualsiasi, sono necessarie tre rotazioni: - una prima rotazione intorno all’asse z, detto angolo di prora o Heading (H), che permette di passare dalla terna 1 o NED alla terna 2; CORSO DI GIS E 3D CITY MODEL 88  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città - una seconda rotazione intorno all’asse trasversale della terna 2, detto angolo di beccheggio o Pitch (P) che permette di passare dalla terna 2 alla terna 3; - una terza rotazione intorno all’asse longitudinale della terna 3, detto angolo di rollio o Roll (R), che permette di passare dalla terna 3 alla terna 4 o body. x 3≡4 p’ 2 𝑅̇ P 1 EST NORD H 1 H 𝑃̇ 2≡3 R 4 y q’ r’ R 𝐻̇ P 4 z 2 1≡2 DOWN Roll (R) Pitch (P) Heading (H) xb xb zb yb yb Figura 4.47 – Rappresentazione degli angoli di Eulero Come già accennato in precedenza, gli angoli di Eulero descrivono, nel caso specifico, l’attitude del frame body rispetto al frame locale di navigazione, la rotazione roll (R) è nota come bank, la rotazione pitch (P) come elevation e, infine la rotazione heading (H) come azimut. Alcuni autori utilizzano il termine attitude per descrivere solo il bank e l’elevation, escludendo l’heading; il bank e CORSO DI GIS E 3D CITY MODEL 89  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città l’elevation sono conosciuti come tilt (Groves, 2008). In questo lavoro per attitude si intenderanno sempre tutte e tre i parametri di assetto. Dal punto di vista analitico, le suddette rotazioni sono rappresentate dalla matrice dei coseni direttori 𝑪𝑏𝑎 = 𝑪𝑏𝑛′ ossia: 𝑪𝑏𝑛′ = 𝑪𝑥 (𝑅)𝑪𝑦 (𝑃)𝑪𝑧 (𝐻 ) [4.47] dove 𝑛′ indica la terna di navigazione NED. Inserendo nella [3.43] le espressioni degli angoli di attitude, si ricava: 𝑐𝑜𝑠𝑃𝑐𝑜𝑠𝐻 𝑐𝑜𝑠𝑃𝑠𝑖𝑛𝐻 −𝑠𝑖𝑛𝑃 𝑪𝑏𝑛′ = (𝑠𝑖𝑛𝑅𝑠𝑖𝑛𝑃𝑐𝑜𝑠𝐻 − 𝑐𝑜𝑠𝑅𝑠𝑖𝑛𝐻 𝑠𝑖𝑛𝑅𝑠𝑖𝑛𝑃𝑠𝑖𝑛𝐻 + 𝑐𝑜𝑠𝑅𝑐𝑜𝑠𝐻 𝑠𝑖𝑛𝑅𝑐𝑜𝑠𝑃 ) [4.48] 𝑐𝑜𝑠𝑅𝑠𝑖𝑛𝑃𝑐𝑜𝑠𝐻 + 𝑠𝑖𝑛𝑅𝑠𝑖𝑛𝐻 𝑐𝑜𝑠𝑅𝑠𝑖𝑛𝑃𝑠𝑖𝑛𝐻 − 𝑠𝑖𝑛𝑅𝑐𝑜𝑠𝐻 𝑐𝑜𝑠𝑅𝑐𝑜𝑠𝑃 che in termini matriciali: 𝐶11 𝐶12 𝐶13 𝑪𝑏𝑛′ = (𝐶21 𝐶22 𝐶23 ) [4.49] 𝐶31 𝐶32 𝐶33 dove ciascuno dei nove componenti rappresenta un coseno direttore, ovvero il coseno che l’asse i della terna body forma con l’asse j nella terna NED. La matrice dei coseni direttori che mette in relazione due terne aventi ciascuna un proprio orientamento, rimane costante se le due terne restano fisse l’una rispetto all’alta. Viceversa, se una terna è in movimento, come nel caso aereo, la matrice varierà nel tempo in funzione della velocità angolare con cui avviene il moto. Pertanto, considerata 𝜔 la velocità angolare con cui la terna (b) si muove rispetto alla terna (a) in un intervallo di tempo infinitesimo 𝑑𝑡, si può ritenere che la rotazione 𝑑𝜃 avvenga intorno ad un asse istantaneo di rotazione tale che: 𝑑𝜃 = 𝜔𝑑𝑡 [4.50] Se 𝑪𝑏𝑎 (𝑡) è la matrice dei coseni direttori relativa ad un istante 𝑡, per un istante successivo 𝑡 + 𝑑𝑡 si ha: 𝑪𝑏𝑎 (𝑡 + 𝑑𝑡 ) = 𝑪 𝑪𝑏𝑎 (𝑡) [4.51] dove C è una matrice che, a sua volta, è il prodotto di tre matrici ciascuna relativa alle tre componenti infinitesime 𝑑𝜃1 , 𝑑𝜃2 , 𝑑𝜃3 del vettore 𝑑𝜃; considerato che tali valori sono piccoli, può porsi: 𝑐𝑜𝑠𝑑𝜃𝑖 ≅ 1 𝑠𝑖𝑛𝑑𝜃𝑖 ≅ 𝑑𝜃𝑖 𝑠𝑖𝑛𝑑𝜃𝑖 𝑠𝑖𝑛𝑑𝜃𝑗 ≅ 0 i quali sostituti nelle equazione [3.43], consentono di ottenere la matrice 𝐶: CORSO DI GIS E 3D CITY MODEL 90  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città 1 𝑑𝜃3 0 1 0 −𝑑𝜃2 1 0 0 𝑪 = (−𝑑𝜃3 1 0) ( 0 1 0 ) (0 1 𝑑𝜃1 ) [4.52] 0 0 1 𝑑𝜃2 0 1 0 −𝑑𝜃1 1 che, in forma matriciale, assume la seguente espressione: 𝑪 = 𝑰 − 𝑨(𝑑𝜃) [4.53] Sostituendo 𝑪 nella [3.52], dividendo per 𝑑𝑡 e tenendo presente la definizione di derivata si ha (Nastro, 2004): 𝑪̇𝑏𝑎 (𝑡 ) = −𝑨(𝜔) 𝑪𝑏𝑎 (𝑡) [4.54] Effettuando la trasposta del primo e del secondo membro, l’equazione [4.54] diviene: 𝑪̇𝑎𝑏 (𝑡) = 𝑪𝑎𝑏 (𝑡)𝑨(𝜔) [4.55] Ricordando che 𝑪𝑖𝑏 è la matrice dei coseni direttori che lega la terna body con quella ECI, la precedente relazione consente di ricavare: 𝑪̇𝑖𝑏 = 𝑪𝑖𝑏 𝑨(𝜔1 ) [4.55] 𝑏 dove 𝜔1 = 𝜔𝑖𝑏 rappresenta la velocità angolare con cui la terna body ruota rispetto a quella inerziale. ′ La matrice 𝐶𝑖𝑛 invece, che lega la terna NED con quella ECI, è: ′ ′ 𝑪̇𝑛𝑖 = −𝑨(𝜔2 )𝑪𝑛𝑖 [4.56] 𝑏 ove 𝜔2 = 𝜔𝑖𝑛′ rappresenta la velocità angolare con cui la terna NED ruota rispetto a quella inerziale. ′ La matrice 𝐶𝑏𝑛 può essere scritta come il prodotto di due matrici, ovvero: ′ ′ 𝑪𝑛𝑏 = 𝑪𝑛𝑖 𝑪𝑖𝑏 [4.57] da cui effettuando la derivata e tenendo conto delle relazioni [4.55] e [4.56], l’equazione precedente diviene: ′ ′ ′ 𝑪̇𝑛𝑏 = 𝑪𝑛𝑏 𝑨(𝜔1 ) − 𝑨(𝜔2 )𝑪𝑛𝑏 [4.58] Si ottengono così nove equazioni differenziali, che risolte danno le componenti della DCM (Nastro, 2004): 𝑡 [4.59] ̇ 𝑑𝜏 𝐶𝑖𝑗 = 𝐶𝑖𝑗 (0) + ∫ 𝐶𝑖𝑗 0 Infine, dall’uguaglianza delle matrici [3.44] e [3.45] si ottiene la relazione che lega gli angoli (RPH) alla matrice dei coseni direttori, ovvero: CORSO DI GIS E 3D CITY MODEL 91  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città 𝐶32 𝑅 = 𝑎𝑟𝑐𝑡𝑎𝑛 ( ) 𝐶33 −𝐶31 𝑃 = 𝑎𝑟𝑐𝑡𝑎𝑛 [ ] [4.60] 2 √1−𝐶31 𝐶21 𝐻 = 𝑎𝑟𝑐𝑡𝑎𝑛 ( ) { 𝐶11 In aerofotogrammetria le tre rotazioni vengono indicate con omega, phi, e kappa (𝜔 , 𝜙, 𝑘 ) e la relazione che li lega con roll, pitch e heading (𝑅, 𝑃, 𝐻 ) è (Dwaik , 1998): 𝜔 = −𝑅 𝜙 = 𝑃 { 𝜋 [4.61] 𝑘 = −𝐻 2 Ottenuta la DCM si conosce, istante per istante, come varia la terna body rispetto a quella di navigazione e, di conseguenza, è possibile ricavare l’evoluzione del moto nel tempo; pertanto risultano derivabili i parametri necessari alla meccanizzazione del sistema inerziale. 4.4.7. MECCANIZZAZIONE DEL SISTEMA INERZIALE Il moto di un mobile nel sistema di riferimento ENU (El-Sheimy, 2004) è modellato dalle seguenti equazioni: 𝑟̅̇ 𝑛 𝑫 𝑣̅ 𝑛 ( 𝑣̅̇ 𝑛 ) = (𝐶𝑏 𝑓 ̅ − (2𝛀𝑖𝑒 + 𝛀𝑒𝑛 )𝑣̅ + 𝑔̅ ) 𝑛 𝑏 𝑛 𝑛 𝑛 𝑛 [4.62] 𝐶𝑏̇ 𝑛 𝑪𝑛𝑏 (Ω𝑏𝑖𝑏 − Ω𝑏𝑖𝑛 ) la prima delle quali rappresenta la relazione tra le coordinate geografiche e le componenti della velocità: 𝑟̅̇ 𝑛 = 𝑫 𝑣̅ 𝑛 dove: 𝑟̅ 𝑛 posizione di un mobile nel frame ENU, espresso in coordinate geografiche; 𝑣̅ 𝑛 velocità nel frame ENU, scomposta nelle tre direzioni East, North ed Up; 𝐷 matrice di trasformazione delle componenti della velocità da lineare ad angolare: 1 0 0 𝜚𝑀 + ℎ 𝑫= 1 0 0 (𝑁𝑃𝑉 + ℎ)𝑐𝑜𝑠𝜑 ( 0 0 1) CORSO DI GIS E 3D CITY MODEL 92  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città avendo indicato con 𝜚𝑀 il raggio di curvatura del meridiano e con 𝑁𝑃𝑉 il raggio di curvatura del primo verticale. La seconda equazione rappresenta, invece, la meccanizzazione della velocità nel frame ENU: 𝑣̅̇ 𝑛 = 𝑪𝑛𝑏 𝑓 ̅𝑏 − (2𝛀𝑛𝑖𝑒 + 𝛀𝑛𝑒𝑛 )𝑣̅ 𝑛 + 𝑔̅ 𝑛 con: 𝑣̅̇ 𝑛 vettore accelerazione nel frame ENU; 𝛀𝑛𝑒𝑛 𝑣̅ 𝑛 accelerazione centripeta relativa al moto nel frame ENU rispetto al frame ECEF; 2𝛀𝑛𝑖𝑒 𝑣̅ 𝑛 accelerazione di Coriolis; 𝑔̅ 𝑛 vettore gravità; 𝑓 ̅𝑏 somma di tutte le accelerazioni nel sistema di riferimento body; 𝜴𝑛𝑒𝑛 𝑛 matrice antisimmetrica del vettore 𝜔𝑒𝑛 (velocità angolare nel sistema di riferimento ENU rispetto al sistema di riferimento ECEF): −𝑉𝑛 𝑉𝑛 𝑉𝐸 𝑡𝑎𝑛𝜑 𝑇 ω𝑛𝑒𝑛 =[ ] 𝜚𝑀 + ℎ 𝑁𝑃𝑉 + ℎ 𝑁𝑃𝑉 + ℎ 𝜴𝑛𝑖𝑒 matrice antisimmetrica del vettore ω𝑛𝑖𝑛 , velocità angolare della Terra: 𝑛 𝜔𝑒𝑛 = [0 𝜎𝑐𝑜𝑠𝜑 𝜎𝑠𝑖𝑛𝜑] La terza equazione rappresenta la meccanizzazione del sistema inerziale: 𝐶𝑏̇ 𝑛 = 𝐶𝑏𝑛 (Ω𝑏𝑖𝑏 − Ω𝑏𝑖𝑛 ) dove 𝛀bib è la matrice antisimmetrica del vettore ω𝑏𝑖𝑏 (velocità angolare del frame body rispetto al frame ECI), nel sistema di riferimento body ed è misurata dai giroscopi, mentre 𝛀𝑏𝑖𝑛 è la matrice antisimmetrica del vettore ω𝑏𝑖𝑛 che rappresenta la velocità angolare del sistema di riferimento di navigazione (NED), rispetto al frame inerziale, nel sistema di riferimento body. Come tutte le misure reali, anche quelle fornite dal sistema inerziale sono affette da errori, per tale motivo si analizzano, nel seguito, le perturbazioni delle diverse variabili che concorrono alla determinazione della posizione e assetto di un mobile. 4.4.8. ERRORI DEL SISTEMA INERZIALE L’insieme delle variabili, posizione, velocità ed assetto ad un generico istante (𝑡), viene definito “stato del sistema inerziale”, indicato con 𝑋(𝑡) e descritto mediante nove componenti: - latitudine; CORSO DI GIS E 3D CITY MODEL 93  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città - longitudine; - quota; - velocità lungo Est; - velocità lungo Nord; - velocità verso l’alto; - rotazione intorno all’asse E-W; - rotazione intorno all’asse N-S; - rotazione intorno all’asse verticale. Lo stato del sistema varia nel tempo e, indicando con 𝑈(𝑡) il vettore delle perturbazioni del sistema, tale variazione può essere descritta da un’equazione differenziale del tipo (Nastro, 2004): 𝑋̇ (𝑡 ) = 𝑓 [𝑋(𝑡 ), 𝑈(𝑡)] [4.63] Purtroppo le componenti non sono note istante per istante, in quanto affette da errori dei quali non si conosce l’andamento. Pertanto si rende necessario ipotizzare un modello matematico che si avvicina a quello reale, detto “stato nominale” 𝑋 ∗ (𝑡), la cuivariazione nel tempo è: 𝑋 ̇∗ (𝑡) = 𝑓 [𝑋 ∗ (𝑡), 𝑈 ∗ (𝑡)] [4.64] con 𝑈 ∗ (𝑡) vettore nominale delle variabili di ingresso. Arrestando ai termini del primo ordine lo sviluppo in serie di Taylor della precedente relazione, si ricava la variazione dello stato del sistema in funzione delle condizioni nominali: 𝜕𝑓 𝜕𝑓 [4.65] 𝑋̇ (𝑡 ) = 𝑋̇ ∗ (𝑡) + ( ) 𝛿𝑋 + ( ) 𝛿𝑈 𝜕𝑋 𝑋=𝑋 ∗ 𝜕𝑈 𝑈=𝑈 ∗ Ponendo: 𝑋(𝑡 ) − 𝑋̇ ∗ (𝑡 ) = 𝛿𝑋(𝑡) = 𝑥 (𝑡 ) 𝑈(𝑡 ) − 𝑈̇ ∗ (𝑡) = 𝛿𝑈(𝑡) = 𝑢(𝑡 ) 𝜕𝑓 ( ) = 𝐹(𝑡) [4.66] 𝜕𝑋 𝑋=𝑋 ∗ 𝜕𝑓 ( ) = 𝐺(𝑡) { 𝜕𝑈 𝑈=𝑈 ∗ e sostituendo le espressioni [3.66] nello sviluppo in serie si ha: 𝑥̇ (𝑡 ) = 𝑭(𝑡 ) 𝑥(𝑡 ) + 𝑮(𝑡 )𝑢(𝑡) [4.67] con: 𝑥 (𝑡 ) vettore le cui componenti rappresentano le deviazioni rispetto allo stato nominale; 𝐹 (𝑡 ) matrice 𝑛 × 𝑛 che esprime la dinamica di evoluzione degli errori; CORSO DI GIS E 3D CITY MODEL 94  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città 𝐺 (𝑡 )𝑢(𝑡) errori generati internamente al sistema inerziale (derive dei giroscopi, errori degli accelerometri). Riscrivendo l’equazione [3.67], rispetto all’equazione della meccanizzazione [3.58], si ottiene la seguente espressione (El-Sheimy, 2004): 𝛿𝑟̅̇ 𝑛 𝐹𝑟𝑟 𝐹𝑟𝑣 𝐹𝑟𝜀 𝛿𝑟̅ 𝑛 03𝑥1 [𝛿𝑣̅̇ ] = [𝐹𝑣𝑟 𝑛 𝐹𝑣𝑣 𝐹𝑣𝜀 ] [𝛿𝑣̅ 𝑛 ] + [ 𝑪𝑛𝑏 𝛿𝑓 ̅𝑏 ] [4.68] 𝜀̅̇𝑛 𝐹𝜀𝑟 𝐹𝜀𝑣 𝐹𝜀𝜀 𝛿𝐶 𝑛 𝑪𝑛𝑏 𝛿𝜔̅𝑏 in cui: 𝜀̅̇𝑛 vettore dei parametri di assetto nel frame ENU; 𝛿𝑓 ̅𝑏 vettore degli errori legati agli accelerometri; ̅𝑏 𝛿𝜔 vettore degli errori nella misura della velocità angolare; −𝑉𝑁 0 0 (𝜚𝑀 + ℎ)2 𝐹𝑟𝑟 = 𝑉𝐸 𝑠𝑖𝑛𝜑 −𝑉𝐸 0 (𝑁𝑃𝑉 + ℎ)𝑐𝑜𝑠 2 𝜑 (𝑁𝑃𝑉 + ℎ)𝑐𝑜𝑠𝜑 ( 0 0 0 ) 1 0 0 𝜚𝑀 + ℎ 𝐹𝑟𝑣 = 1 0 0 (𝑁𝑃𝑉 + ℎ)𝑐𝑜𝑠𝜑 ( 0 0 1) 𝐹𝑟𝜀 = 03𝑥3 𝑉𝐸 𝑉𝑁 𝑉𝐸 𝑉𝑁 𝑉𝐸 𝑉𝑁 𝑡𝑎𝑛𝜑 2𝜎(𝑉𝑈 𝑠𝑖𝑛𝜑 + 𝑉𝑁 𝑐𝑜𝑠𝜑) + 0 − (𝑁𝑃𝑉 + ℎ)𝑐𝑜𝑠 2 𝜑 (𝑁𝑃𝑉 + ℎ) 2 (𝑁𝑃𝑉 + ℎ)2 𝑉𝐸 2 𝑉𝑁 𝑉𝑈 𝑉𝐸 2 𝑡𝑎𝑛𝜑 𝐹𝑣𝑟 = −2𝜎 𝑉𝐸 𝑐𝑜𝑠𝜑 − 0 − (𝑁𝑃𝑉 + ℎ)𝑐𝑜𝑠 2 𝜑 (𝜚𝑀 + ℎ)2 (𝑁𝑃𝑉 + ℎ)2 𝑉𝐸 2 𝑉𝑁 2 2𝑔 −2𝜎𝑉𝐸 𝑠𝑖𝑛𝜑 0 − − + ( (𝑁𝑃𝑉 + ℎ) 2 (𝜚𝑀 + ℎ) 2 𝑅 + ℎ) 𝑉𝑁 𝑡𝑎𝑛𝜑 𝑉𝑈 𝑉𝐸 𝑡𝑎𝑛𝜑 𝑉𝐸 − 2𝜎 𝑠𝑖𝑛𝜑 + 2𝜎 𝑐𝑜𝑠𝜑 − 𝑁𝑃𝑉 + ℎ 𝑁𝑃𝑉 + ℎ 𝑁𝑃𝑉 + ℎ 𝑁𝑃𝑉 + ℎ 2 𝑉𝐸 𝑡𝑎𝑛𝜑 𝑉𝑈 𝑉𝑁 𝐹𝑣𝑣 = −2𝜎 𝑉𝐸 𝑠𝑖𝑛𝜑 − − − 𝑁𝑃𝑉 + ℎ 𝜚𝑀 + ℎ 𝜚𝑀 + ℎ 2 𝑉𝐸 𝑉𝑁 2𝜎 𝑐𝑜𝑠𝜑 + 0 ( 𝑁𝑃𝑉 + ℎ 𝜚𝑀 + ℎ ) 0 𝑓𝑈 −𝑓𝑁 𝐹𝑣𝜀 = (−𝑓𝑈 0 𝑓𝐸 ) 𝑓𝑁 −𝑓𝐸 0 𝑉𝑁 0 0 − (𝜚𝑀 + ℎ)2 𝑉𝐸 𝐹𝜀𝑟 = 𝜎 𝑠𝑖𝑛𝜑 0 (𝑁𝑃𝑉 + ℎ)2 𝑉𝐸 𝑉𝐸 𝑡𝑎𝑛𝜑 −𝜎 𝑐𝑜𝑠𝜑 − 0 ( (𝑁𝑃𝑉 + ℎ)𝑐𝑜𝑠 2 𝜑 (𝑁𝑃𝑉 + ℎ)2 ) CORSO DI GIS E 3D CITY MODEL 95  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città 1 0 0 𝜚𝑀 + ℎ 1 𝐹𝜀𝑣 = − 0 0 𝑁𝑃𝑉 + ℎ 𝑡𝑎𝑛𝜑 − 0 0 ( 𝑁𝑃𝑉 + ℎ ) 𝑉𝐸 𝑡𝑎𝑛𝜑 𝑉𝐸 0 𝜎 𝑠𝑖𝑛𝜑 + −𝜎 𝑐𝑜𝑠𝜑 − 𝑁𝑃𝑉 + ℎ 𝑁𝑃𝑉 + ℎ 𝑉𝐸 𝑡𝑎𝑛𝜑 𝑉𝑁 𝐹𝜀𝜀 = −𝜎 𝑠𝑖𝑛𝜑 − 0 − 𝑁𝑃𝑉 + ℎ 𝜚𝑀 + ℎ 𝑉𝐸 𝑉𝑁 𝜎 𝑐𝑜𝑠𝜑 + 0 ( 𝑁𝑃𝑉 + ℎ 𝜚𝑀 + ℎ ) Dai parametri di output del sistema inerziale, e dai relativi errori, si è in grado di stabilire posizione e traiettoria di un mobile; tuttavia l’entità degli errori può risultare non compatibile con le precisioni attese nelle applicazioni geomatiche a causa della deriva nel tempo dei giroscopi. Per tale motivo, è necessario integrare le informazioni ottenute dai sistemi inerziali con ulteriori dati provenienti da fonti esterne come, ad esempio, quelli da sistema GNSS. 4.4.9. LA NAVIGAZIONE INTEGRATA Allo scopo di migliorare le informazioni fornite dal sistema inerziale si ricorre alla navigazione integrata, cioè ad un sistema idoneo al trattamento di una serie di misure, di diversa provenienza, tutte concorrenti alla determinazione ottimale della posizione e dei parametri di guida di un aereo (Nastro, 2004). Si dimostra, infatti, come due serie di misure omogenee e tra loro indipendenti, ottenute con metodi diversi, consentono di ricavare la stima ottimale di una quantità scalare 𝑥, ovvero individuando idonei strumenti matematici. Tra i numerosi strumenti matematici che possono essere utilizzati per la stima stocastica di misure di sensori rumorosi, uno dei più noti, e di uso più frequente, è il Filtro di Kalman (FK) così detto dopo che, nel 1960, Rudolph E. Kalman pubblicò il suo famoso documento descrivente una soluzione ricorsiva del problema del filtraggio lineare discreto dei dati (Kalman, 1960); dal momento della sua introduzione, esso è stato oggetto di approfondite ricerche e applicazioni. Il Filtro di Kalman è un insieme di equazioni matematiche che implementano uno stimatore predittore-correttore che è ottimale poiché minimizza, quando sono soddisfatte alcune presunte condizioni, la covarianza dell’errore stimato. I cinque elementi fondamentali dell’ FK (Figura 3.13) sono: il vettore di stato e la covarianza, il modello di sistema, il vettore delle misurazione e la covarianza, il modello delle misurazione e l'algoritmo; nella figura 4.25 la linea continua indica CORSO DI GIS E 3D CITY MODEL 96  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città i flussi di dati che sono sempre presenti, mentre quella tratteggiata i quelli presenti solo in alcune applicazioni (Groves, 2008). Figure 4.25 - Elementi del Filtro di Kalman In virtù di tali elementi, da una prima serie di stime, il FK opera in modo ricorsivo aggiornando le proprie stime come media ponderata dei valori precedenti e dei nuovi valori derivati dagli ultimi dati di misura, ciò differentemente dagli algoritmi di stima non ricorsivi che derivano le stime dei parametri da tutto l'insieme di dati di misura senza stime a priori. Per le applicazioni in tempo reale, come la navigazione, l'approccio ricorsivo è un processore più efficiente poiché, a ogni iterazione, solo i nuovi dati di misura devono essere necessariamente trattati, mentre i vecchi possono essere scartati; per tale motivo il FK trova particolare applicazione nel settore della navigazione autonoma o assistita e viene adottato nell’integrazione GNSS/INS. Esistono diverse “architetture” di integrazione GPS/INS che possono essere classificate come descritto in seguito (Petovello 2003; Cazzaniga, 2007). Nell’integrazione uncoupled le stime GPS sono usate semplicemente per “re-inizializzare” posizione e velocità dell’INS a intervalli regolari di tempo. Nell’integrazione loosely coupled (LC) le elaborazioni dei dati sono disposte in cascata: prima si ha un filtro GPS che fornisce posizione e velocità, successivamente si valutano le differenze tra la stima GPS e quella INS e si stimano gli errori di quest’ultima usati per correggere la soluzione inerziale tramite le equazioni variazionali, il cui schema dell’algoritmo viene riportato in figura 3.15. Il metodo di integrazione loosely contiene due filtri (uno per le misure GPS e uno per quelle INS) che operano in modo indipendente l’uno dall’altro: 1. dal processamento, con un filtro di “Kalman GPS”, delle misure GPS raw si n n ( determinano la posizione e la velocità da satellite rGPS , vGPS ; ) 2. processando, con le equazioni di meccanizzazione, le misure INS raw  ibb , v b  si CORSO DI GIS E 3D CITY MODEL 97  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città determinano la posizione e la velocità da INS  rINS n n , vINS ; 3. la differenza r n , v n  tra posizione e velocità GPS rGPS n n , vGPS  e INS  rINSn , vINS n  (ottenute in 1. e 2.) viene data in input ad un filtro “Kalman INS”; 4. tale filtro stima gli errori di posizione e velocità  r n ,  vn  nonché l’errore di disallineamento   n  ; 5. gli errori stimati vengono impiegati per aggiornare la posizione e la velocità r n INS n , vINS  (ricavata in 2.) al fine di ottenere un vettore di stato completo r , v , R n n n b   Cbn. Figura 4.26 – Diagramma di flusso dell’integrazione loosely Il vantaggio del metodo d’integrazione loosely è la sua semplicità nella realizzazione nonché la sua robustezza; infatti se uno dei sensori (INS o GPS) fallisce, la soluzione sarà data dall’altro sensore). Un altro vantaggio è la rapidità di processamento dell’algoritmo dovuto, generalmente, a vettori di stato più piccoli (Petrovello, 2003). Tra gli svantaggi va invece sottolineato che con una copertura povera (meno di quattro satelliti) è impossibile fornire un aggiornamento della misura dal filtro GPS. Nella strategia di integrazione tightly coupled (TC) le misure di “pseudorange” e “pseudorange rate” fornite dal GPS, sono confrontate con le stime di queste quantità generate a partire dalle misure inerziali. L’algoritmo tightly coupled, contiene un solo un filtro di Kalman (Filtro di Kalman GPS/INS) e può essere schematizzato (Figura 3.16) nel seguente modo: 1. dal processamento, con le equazioni di meccanizzazione, delle misure INS raw  b ib  , v b si ricava la posizione e la velocità da INS rINS n  n , vINS ;  CORSO DI GIS E 3D CITY MODEL 98  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città 2. l’unione di tali valori  rINS n n , vINS  con le effemeridi GPS raw, consente di predire le misure di pseudorange e Doppler INS ,INS ; 3. parallelamente, processando i dati GPS raw, si ottengono misure di range e Doppler GPS GPS ,GPS ; 4. la differenza  ,  tra le misure di pseudorange e Doppler predette INS ,INS  e quelle misurate GPS ,GPS  (ottenute in 2. e 3.) viene data in input ad un filtro di Kalman GPS/INS; 5. tale filtro stima gli errori di posizione e velocità  r n ,  vn  nonché l’errore di disallineamento  n  ; 6. gli errori stimati vengono impiegati per aggiornare la posizione e la velocità  rINS n n , vINS  (ricavata in 1.) al fine di ottenere un vettore di stato completo  r n , vn , Rbn  Cbn . Figura 4.27 – Diagramma di flusso dell’integrazione tightly A causa dell’utilizzo dei due differenti filtri, tale architettura viene anche denominata “metodo di integrazione decentralizzato”. Infine, nella deeply coupled o ultra-tightly coupled si combina il tracciamento del segnale GPS e l’integrazione GPS/INS in un unico filtro di Kalman che viene impiegato per analizzare i dati di posizionamento sia all’indietro “backword” che in predizione mediante un algoritmo di smoothing. CORSO DI GIS E 3D CITY MODEL 99  CAPITOLO IV – Sensori passivi e attivi aviotrasportati per la rappresentazione 3D delle città 4.4.10. La navigazione integrata applicata nel campo aerofotogrammetrico Un sistema di navigazione integrata (Positioning and attutude – POS) connessa ad una fotocamera permette di conoscere in ogni momento la posizione e l'assetto del sensore grazie all'integrazione di sistemi GNSS e INS. I vantaggi di un'integrazione GNSS/INS sono che le stime INS possono essere corrette dai dati GNSS e che l'INS può fornire aggiornamenti della posizione e d