Modelado y simulación de un receptor digital Loran

Equation Chapter 1 Section 1 Trabajo Fin de Grado Grado en Ingeniería Aeroespacial Int. en Navegación Aérea Modelado y Simulación de un Receptor Digital Loran Autor: Eduardo de la Rosa Padrón Tutor: Joaquín Granado Romero Dep. Ingeniería Electrónica Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2017 Trabajo Fin de Grado Grado en Ingeniería Aeroespacial Int. en Navegación Aérea Modelado y Simulación de un Receptor Digital Loran Autor: Eduardo de la Rosa Padrón Tutor: Joaquín Granado Romero Profesor Titular de Universidad Dep. Ingeniería Electrónica Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2017 Trabajo Fin de Grado: Modelado y Simulación de un Receptor Digital Loran Autor: Eduardo de la Rosa Padrón Tutor: Joaquín Granado Romero El tribunal nombrado para juzgar el Proyecto arriba indicado, compuesto por los siguientes miembros: Presidente: Vocales: Secretario: Acuerdan otorgarle la calificación de: Sevilla, 2017 El Secretario del Tribunal A los que creyeron en mí cuando yo no lo hacía: mi familia Agradecimientos Lo primero es lo primero: la familia, sin los que todo esto no habría sido posible. Gracias al esfuerzo de todos ellos, tanto de los que viven bajo mi mismo techo como de los que no, hoy yo no estaría aquí. Soy quien y como soy porque ellos han tenido la paciencia de enseñarme, al golpito, cómo ser una persona humilde, agradecida y honrada. Me gustaría también expresar lo afortunado que me siento por haber podido realizar este proyecto con aquel profesor con quien siempre deseé hacerlo, Joaquín Granado. Querría agradecerle de forma especial y efusiva el haber aceptado de manera incondicional ser mi tutor, por toda su ayuda prestada, su trato cercano y cómodo, su consejo, y por su dedicación. Dicen que los amigos son la familia que uno elige. Yo he elegido una familia más bien reducida, y que ahora anda un tanto dispersa, pero con la seguridad de que hemos establecido lazos no volátiles. Debo pedirles perdón y darles las gracias por el mismo motivo: soportarme durante estos años compartidos donde cada uno ha aportado su hilo particular para confeccionar un saco de recuerdos inmenso. De entre ellos debo nombrar a Mauricio Di Rado por haberme ayudado de forma completamente desinteresada con el ámbito idiomático de este trabajo. No puedo finalizar estos agradecimientos sin hacer una mención especial a Ruymán Cruz, profesor del Departamento de Análisis Matemático de la Universidad de La Laguna, pero ante todo amigo, por su valioso apoyo a lo largo de todos estos años, tanto en lo personal como en lo académico. Sus meticulosas aportaciones han contribuido a aportar el rigor que este documento se merece. Eduardo de la Rosa Padrón Sevilla, 2017 i 3.4. Estimación de la posición . . . . . . . . . . . . . . . . . . . . . . . . . 41 3.4.1. Fuentes de error . . . . . . . . . . . . . . . . . . . . . . . . . . 43 3.4.2. Fuentes adicionales . . . . . . . . . . . . . . . . . . . . . . . . 47 4. Resultados de las simulaciones 49 4.1. Metodología................................ 49 4.2. Simulaciones realizadas . . . . . . . . . . . . . . . . . . . . . . . . . . 50 4.2.1. Posicionamiento estático . . . . . . . . . . . . . . . . . . . . . 50 4.2.2. Posicionamiento móvil . . . . . . . . . . . . . . . . . . . . . . 55 4.3. Efectodelcanal.............................. 66 4.3.1. Ruido ............................... 66 4.3.2. Atenuación ............................ 69 5. Conclusiones y trabajos futuros 73 5.1. Conclusiones................................ 73 5.2. Propuestas de ampliación . . . . . . . . . . . . . . . . . . . . . . . . 74 Nomenclatura 77 Bibliografía 81 Anexo A. Información de las Loran-C chain 85 Anexo B. Método de Chan 91 B.1. Modelo matemático general para la estimación hiperbólica de la posición 91 B.2. Propuesta de resolución del Método de Chan . . . . . . . . . . . . . . 92 viii Índice de figuras 2.1. El sistema Loran-C: dos transmisores cada uno de ellos radiando breves pulsos de anchura especificada con un intervalo de repetición fijo . 6 2.2. Línea LOP producida por dos estaciones transmisoras emitiendo pulsossimultáneamente ........................... 6 2.3. Familia de LOPs producidas por dos estaciones transmisoras emitiendo pulsos simultáneamente . . . . . . . . . . . . . . . . . . . . . . . . 7 2.4. Modificación de las LOPs de la Fig. 2.3. La estación B no comienza su transmisión hasta que es disparada por un pulso proveniente de la estaciónA................................. 8 2.5. Una modificación adicional en las LOPs de la Fig. 2.3. La estación B no sólo debe esperar la recepción de un pulso de la estación A sino además un CD, alterando el valor de la diferencia temporal que caracterizaacadaLOP.......................... 9 2.6. Cronograma de emisión y recepción de la señal . . . . . . . . . . . . . 9 2.7. Un chain está constituido por un master y hasta un total de cinco slaves ................................... 11 2.8. Forma del pulso unitario del sistema Loran-C . . . . . . . . . . . . . 13 2.9. Transmisión Loran-C por single pulsing ................. 14 2.10. Transmisión Loran-C por multiple pulsing ............... 14 2.11. Restricciones para la asignación del ED . . . . . . . . . . . . . . . . . 15 2.12. Ground waves ysky waves ........................ 19 2.13. Área típica de cobertura para el día y la noche del sistema Loran-C . 20 2.14. Diagrama funcional del equipamiento de una estación Loran-C . . . . 21 2.15. Disposición típica de los equipos de transmisión . . . . . . . . . . . . 22 2.16. Ejemplos de receptores Loran-C embarcados en aeronaves . . . . . . . 24 2.17. NGA Loran chart para el chain 9940, parejas 9940-X y 9940-Y . . . . 24 3.1. Diagrama de la arquitectura del receptor . . . . . . . . . . . . . . . . 25 3.2. Forma directa de un filtro FIR discreto de orden N.......... 26 3.3. Módulo de la respuesta en frecuencia de los filtros FIR BPF y LPF . 27 3.4. PSD de las señales entrante y saliente a los filtros BPF y LPF, así como el de su original transmitida . . . . . . . . . . . . . . . . . . . . 28 3.5. CircuitoDdEdigital ........................... 29 3.6. Diferentes comportamientos del DdE ante distintos valores de τ1yτ231 3.7. Comportamiento del DdE para los parámetros τ1yτ2ajustados . . . 32 3.8. Umbrales, medidas del ruido e instantes CheckMax . . . . . . . . . . 33 ix 3.9. Porcentaje de los valores aleatorios que caen dentro de la banda centrada en la media de una distribución normal y con un ancho de dos, cuatro y seis desviaciones estándar . . . . . . . . . . . . . . . . . . . 35 3.10. Posibles escenarios que pueden presentarse a la hora de validar un máximo en el instante CheckMax . . . . . . . . . . . . . . . . . . . . 36 3.11. OAs antes y después del procedimiento de retraso . . . . . . . . . . . 38 3.12. Casuística presentable en la detección de estaciones . . . . . . . . . . 39 3.13. Diagrama de flujo detallado del control del receptor . . . . . . . . . . 40 3.14. Señales producidas de cada etapa del receptor . . . . . . . . . . . . . 42 3.15. Error de estimación de la posición debida al movimiento del vehículo 43 3.16. Señal formada en el receptor para diferentes valores de Ts....... 44 3.17. Áreas de incertidumbre en función del gradiente hiperbólico . . . . . 45 3.18. Familia de LOPs de una pareja master-slave donde se ilustran los excesivos gradientes que se tienen cerca de la baseline extension . . . 46 4.1. Esquema metodológico propuesto para el estudio del sistema Loran-C 49 4.2. Comparación entre la (in)exactitud y la (im)precisión . . . . . . . . . 52 4.3. Error DRMS en función del SNR y el gradiente para el posicionamientoestático .............................. 53 4.4. Posiciones estimadas para diferentes valores de SNR durante posicionamientoinmóvil ............................. 55 4.5. Definición gráfica de ATE y XTE . . . . . . . . . . . . . . . . . . . . 56 4.6. Definición de precisión y exactitud en términos probabilísticos . . . . 57 4.7. Posiciones estimadas para diferentes valores de SNR volando a 250 km/h en trayectoria recta . . . . . . . . . . . . . . . . . . . . . . . . 59 4.8. Posiciones estimadas para diferentes valores de SNR volando a 900 km/h en trayectoria recta . . . . . . . . . . . . . . . . . . . . . . . . 59 4.9. Posiciones estimadas diezmadas un factor 12 frente a las originales para diferentes valores de SNR volando a 250 km/h en trayectoria recta 60 4.10. Histogramas y ajuste a distribución normal de errores ATE y XTE para diferentes valores de SNR volando a 250 km/h en trayectoria recta; estimaciones originales . . . . . . . . . . . . . . . . . . . . . . . 61 4.11. Histogramas y ajuste a distribución normal de errores ATE y XTE para diferentes valores de SNR volando a 250 km/h en trayectoria recta; estimaciones diezmadas . . . . . . . . . . . . . . . . . . . . . . 62 4.12. Posiciones estimadas diezmadas un factor 12 frente a las originales para diferentes valores de SNR volando a 900 km/h en trayectoria recta 63 4.13. Histogramas y ajuste a distribución normal de errores ATE y XTE para diferentes valores de SNR volando a 900 km/h en trayectoria recta; estimaciones originales . . . . . . . . . . . . . . . . . . . . . . . 64 4.14. Histogramas y ajuste a distribución normal de errores ATE y XTE para diferentes valores de SNR volando a 900 km/h en trayectoria recta; estimaciones diezmadas . . . . . . . . . . . . . . . . . . . . . . 65 4.15. Posiciones estimadas para diferentes valores de SNR volando a 900 km/henviraje .............................. 67 4.16. Efecto del ruido sobre la señal transmitida para diferentes valores de SNR .................................... 69 x 4.17. Efecto de la atenuación de la señal en la detección/validación de un grupodepulsos.............................. 71 xi Índice de tablas 2.1. Valores que adoptan las funciones pcXY ................ 16 2.2. Permitividad y conductividad de diferentes medios de propagación . . 21 4.1. Parámetros de configuración establecidos para las simulaciones . . . . 51 4.2. Información del chain empleado..................... 51 4.3. Valores estimados de la media y desviación estándar para el ajuste gaussiano de los errores en el caso de vuelo a 250 km/h . . . . . . . . 58 4.4. Valores estimados de la media y desviación estándar para el ajuste gaussiano de los errores en el caso de vuelo a 900 km/h . . . . . . . . 58 4.5. Valores estimados de la media y desviación estándar para el ajuste gaussiano de los errores en el caso de vuelo a 900 km/h . . . . . . . . 66 A.1. Información de las Loran-C chain en coordenadas WGS-84 . . . . . . 86 A.2. Información de las Loran-C chain en coordenadas WGS-84 (cont.) . . 87 A.3. Información de las Loran-C chain en coordenadas WGS-84 (cont.) . . 88 A.4. Información de las Loran-C chain en coordenadas WGS-84 (cont.) . . 89 xiii Capítulo 1 Introducción: justificación y objetivos Lucha, mientras los demás están durmiendo. Aprende, mientras el resto está de fiesta. Vive, como ellos sueñan. Katinka Hosszú, nadadora húngara 1.1. Reseña histórica Elsistema Loran (LOng RAnge Navigation) es un sistema de radionavegación hiperbólico que permite al receptor determinar su ubicación gracias a la escucha de señales de radio de baja frecuencia transmitidas por estaciones fijas en tierra. En general, cualquier tipo de navegador se posiciona gracias a la intersección de dos líneas de posición situadas sobre la superficie terrestre. Cada una de estas líneas es el lugar geométrico en el que un determinado parámetro que las caracterice tiene un valor específico. Es conveniente clasificar los sistemas de ayuda a la navegación acorde con la forma que puede adoptar sus líneas de posición [1]: Posicionamiento mediante la determinación de dos direcciones: las líneas son radios (por ejemplo, el sistema ADF). Posicionamiento mediante la determinación de una dirección y una distancia: por coordenadas polares (por ejemplo, el sistema RADAR). Posicionamiento mediante la determinación de dos distancias: las líneas son circunferencias (por ejemplo, el sistema DME). Posicionamiento mediante la determinación de dos diferencias de distancia: las líneas son hipérbolas, como es el caso del sistema Loran. Por suerte o por desgracia, el ingenio del ser humano se agudiza durante las etapas de guerra con el fin de estar en una mejor posición que el enemigo. Un gran número de descubrimientos e inventos se han logrado durante las épocas bélicas, y Loran no escapa a este macabro origen. Durante la Primera Guerra Mundial se ideó un método para localizar las armas enemigas camufladas, que consistía en medir las diferencias entre los tiempos con los que eran escuchados sus disparos por tres puestos de escucha diferentes. Estos puestos se encontraban eléctricamente conectados, por lo 1 que las tres denotaciones recibidas de un mismo disparo efectuado por un arma podían ser registrados en un cronógrafo. El intervalo de tiempo transcurrido desde la llegada de un sonido a uno de los puestos hasta su llegada a otro es una medida de cuánto la distancia del arma al segundo puesto excede la distancia al primero. Esta diferencia es precisamente la constante que define una hipérbola, en cuyos focos se encuentran los dos puestos de escucha. El intervalo temporal entre la recepción del sonido del disparo por el segundo puesto y un tercer puesto adicional localiza el arma sobre una línea de posición hiperbólica que pasa entre los dos puestos y que es cóncava hacia el más cercano. De manera análoga, se tendría otra línea de posición diferente para los dos primeros puestos. Estas dos líneas de posición hiperbólicas se intersectan en el arma, lo que permitía a los soldados conocer su localización y así poder destruirla. Este es el principio en el que se basan los sistemas de navegación denominados hiperbólicos, cambiando las señales sonoras por señales de radio. Gran Bretaña y Estados Unidos fueron los pioneros en su desarrollo desde su aparición hasta llegar al que fue el sistema hiperbólico por excelencia: Loran-C [1]. Este sistema ha estado operativo en los EE.UU. desde los años 70 [2], proporcionando una precisión en el posicionamiento de centenares de metros. Cuando el GPS (Global Positioning System) apareció en la década de 1980, con una precisión de decenas de metros, muchos empezaron a considerar a Loran-C irrelevante. Sin embargo, otros veían en él una fuente esencial de posicionamiento y tiempo en la que confiar en caso de fallo del GPS. En septiembre de 2001, casi simultáneamente a los atentados del 9/11, el gobierno de los EE.UU. publicó su “Volpe Report” [3]. En él se explicaba el grado en el que tanto los EE.UU. como otros países estaban empezando a basar sus infraestructuras críticas en el sistema GPS. El informe manifestaba la vulnerabilidad del GPS (y otros sistemas GNSS -Global Navigation Satellite System-) a las interrupciones del servicio debidas a interferencias intencionadas y no intencionadas. Se reconoció a Loran como una potencial solución para este importante problema. Esto atrajo el interés mundial y trajo consigo un ímpetu en modernizar el sistema Loran en EE.UU. El programa americano de evaluación y modernización de Loran fructificó en una nueva versión que ofrecía notables mejoras en las prestaciones, versión que se denominó Enhanced Loran (eLoran). Este renovado sistema ofrecía una mayor exactitud, integridad y continuidad a la vez que se preservaban los requisitos de disponibilidad del Loran-C tradicional. Estas mejoras se han obtenido principalmente a través de la adición de un data channel (LDC, Loran Data Channel), siendo esta la principal diferencia entre ambas versiones de Loran. Este data channel permite a eLoran cumplir con los requisitos que se imponen para aterrizajes haciendo uso de instrumentos de aproximación de no precisión (NPA, Non-Precision Approach), y de llevar con seguridad a los barcos a puerto (HEA, Harbor Entrance and Approach) en condiciones de poca visibilidad. eLoran se mostró capaz también de proveer referencias de tiempo y frecuencia extremadamente precisas necesarias para los sistemas de telecomunicaciones que emplean comunicaciones de audio e internet. La idea del LDC es enviar con un bajo bitrate mensajes tales como correcciones diferenciales en tiempo real, avisos, e incluso información sobre la integridad de la señal al receptor del usuario por medio de la transmisión Loran. Estas correcciones son ofrecidas por estaciones de referencia que detectan diminutas variaciones en la señal eLoran, de manera similar a cómo operan las estaciones de referencia de los GNSS diferenciales, permitiendo a los receptores corregir estas variaciones. Ya que eLoran emplea transmisores de gran potencia para 2 emitir señales de baja frecuencia (no microvatios y microondas como GNSS), es muy improbable que pueda ser perturbada o interferida por las mismas causas por las mismas causas que afectan a las señales GNSS. Esto implica que pequeños receptores eLoran, incluso incorporados en las mismas unidades GNSS, pueden mitigar el impacto de las interrupciones de los GNSS [4]. La FAA (Federal Aviation Administration) declaró que “[...] la modernización del sistema Loran podría satisfacer los actuales restricciones NPA, HEA y de tiempo/frecuencia en los EE.UU. y podría ser empleado para mitigar los efectos operacionales de una interrupción en el servicio GPS, permitiendo de esta manera a los usuarios conservar los beneficios derivados del uso del GPS” [5]. Esto no impidió que años más tarde el entonces presidente de EE.UU., Barack Obama, cumpliera su promesa de eliminar programas federales, como el del sistema Loran-C, considerados innecesarios e inviables económicamente a causa de la importante inversión en nueva infraestructura y equipamiento del usuario necesaria para implantar eLoran. La USCG (United States Coast Guard) declaró que “Loran-C [...] se había convertido en un sistema anticuado que ya no era requerido por las fuerzas armadas, el sector del transporte o los intereses de seguridad de la nación”, finalizando así las transmisiones de todas las señales Loran-C en terreno americano el pasado día 8/2/2010 [6]. Consecuentemente, el sistema entró en decadencia a nivel mundial. 1.2. Objetivos Este trabajo surge en relación con un concepto propuesto por la Universidad de Delft conocido como Eurofix [7], en el que se trata de combinar sistemas GNSS (como GPS y Galileo) y Loran-C en un único sistema integrado. El primer paso para ello es analizar individualmente cada sistema y sus actuaciones. Se ha optado por estudiar en concreto el sistema Loran-C por ser bastante más desconocido que los sistemas GNSS, sobre todo en el mundo aeronáutico. El propósito principal de este trabajo será el estudio detallado de este sistema desde el punto de vista de la señal. Para ello primero habrá que comprender tanto los fundamentos que rigen el funcionamiento de este sistema hiperbólico cómo las características propias de la modulación Loran-C que posibiliten generar la señal. Conocidos estos aspectos se está en disposición de proponer una arquitectura digital de receptor que permita la detección de la señal Loran-C y la posterior extracción de la información necesaria para poder estimar la posición del receptor. En último lugar, y con el propósito doble de validar el diseño planteado y conocer las actuaciones del propio sistema de navegación, se realizarán una serie de simulaciones en diferentes escenarios y configuraciones. Se ha procurado hacer hincapié en la contrastación de la información mediante la consulta en multitud de referencias de legítima procedencia que enriquezcan y aporten veracidad. 1.2.1. Metodología e hipótesis generales De manera conceptual, la metodología llevada a cabo partirá de la definición arbitraria de una trayectoria real seguida por el vehículo, en este caso una aeronave. El conocimiento de estos puntos permite conformar una señal de el transmisor que 3 la siguiente manera: TD = 2a c±BTT+CD, (2.6) donde 2aes la constante definida en (2.3), cla velocidad de propagación de la señal, BTT el baseline travel time, CD el coding delay, y donde los signos +y− corresponden a los semiplanos master yslave respectivamente. Estos semiplanos hacen referencia a los dos sectores que quedan a ambos lados de la mediatriz de la baseline: aquel que contenga la posición de la estación master será el semiplano master (HM), mientras que el contenga a la estación slave será el semiplano slave (HS). Por otro lado, de la propia ecuación de la hipérbola se sigue que 0≤2a≤d(F1,F2),(2.7) que en forma temporal puede expresarse según 0≤2a c0≤BTT, (2.8) permitiendo eliminar el valor absoluto en (2.6). Se tiene así que a(TD) =    c0 2(BTT +CD −TD)si TD ∈[CD, CD +BTT], c0 2(TD −BTT −CD)si TD ∈(CD +BTT, CD + 2BTT], (2.9) TD(a) =        BTT +CD −2a c0 si dS(~x)≤dM(~x)⇔~x ∈HS, BTT +CD +2a c0 si dS(~x)> dM(~x)⇔~x ∈HM. (2.10) En las relaciones (2.9) y (2.10) los parámetros c0, BTT y CD son constantes: BTT y CD por ser preestablecidos y c0por ser una constante universal. Los dos bloques principales (Tx y Rx) del sistema Loran-C modelado harán uso de estas relaciones. TD i −* )− iv aii −* )− iii ~x (2.11) Tx: en el conformado de la señal se parte de una posición conocida ~x para llegar a TD. Por lo tanto, habría que aplicar los procesos iii yiv, que corresponden a las ecs. (2.4) y (2.10) respectivamente. Rx: para la estimación de la posición se parte de TD para llegar a (una estimación de) ~x. Para ello se atraviesan los procesos iyii, correspondientes a las ecs. (2.9) y (2.12) respectivamente. x2 a2−y2 c2 0−a2= 1 (2.12) La ec. (2.12) no es más que la ecuación canónica de la hipérbola en coordenadas cartesianas xey. Todos los procesos de (2.11) son unívocos salvo el proceso ii, que es trivial: un punto en el plano es atravesado por una única hipérbola (implicando un único valor de a), pero ella misma atraviesa también una infinidad de otros puntos. Por lo tanto, el conocimiento de un valor de posición ~x brinda un único valor de a(proceso iii), mientras que el conocimiento de un valor de aofrece una inmensidad de posibles puntos. Esta es la razón por la que es necesario disponer de dos curvas hiperbólicas (dos LOPs) con el fin de intersectarlas para hallar un único punto. 10 2.1.2. Arquitectura Loran-C El sistema Loran-C se encuentra constituido por un conjunto (chain) que consta de un mínimo de tres a un máximo de seis estaciones de transmisión. Una de ellas siempre se designa como master (M), mientras que el resto se comportarán como estaciones slaves, designadas victor (V), whisky (W), x-ray (X), yankee (Y) y zulu (Z), acorde con el International Phonetic Alphabet. Existe la posibilidad de que una misma estación master sea compartida por dos chains distintos (dual-rated). Figura 2.7: Un chain está constituido por un master y hasta un total de cinco slaves (Fuente: [8]). Todas las estaciones de un mismo chain se encuentran sincronizadas entre sí, por lo que las señales de los slaves poseen precisas relaciones temporales con las transmisiones del master. Detalles de las estaciones que conforman cada chain, así como de los parámetros característicos de cada una de ellas se encuentran recopilados en el Anexo A. Es evidente que el método más sencillo (entendido como que suponga la menor carga de maniobra para el piloto) para ir de un punto a otro es mediante una ruta con rumbo constante. Por este motivo las aerovías del espacio aéreo se han confeccionado tradicionalmente haciendo uso de radioayudas VOR como sus puntos de origen/destino, donde el piloto únicamente tendrá que interceptar un radial de la radiobaliza y mantenerse en él. Si este sistema de guiado punto a punto quisiera emplearse con Loran-C, su geometría hiperbólica entrañaría cierta complejidad a la hora del guiado por no ser curvas de rumbo constante. Pese a ello, Loran-C ofrece una característica adicional importante: es posible dirigirse a cualquier punto geográfico dentro de la cobertura del chain. Las ventajas de este sistema son obvias en comparación con otros basados en radiobalizas, en los que el modo de vuelo se restringe a una trayectoria desde/hacia la radiobaliza [11]. 11 2.2. Características de la señal 2.2.1. Forma del pulso Loran-C Las transmisiones del sistema Loran-C son transmisiones de tipo ‘pulsadas’. Pulsos de radiofrecuencia de muy corta duración son radiados a intervalos periódicos que son mucho más largos que la duración del pulso. Durante el intervalo entre pulsos sucesivos de una estación existe un período de ‘silencio’ [12]. La forma de estos pulsos se encuentra definida en el documento COMDTINST M16562.4A de la USCG [13]. En primer lugar se define la envolvente de un pulso individual de la forma: E(t) =      0para t < τ A(t−τ)2exp[−2(t−τ)/65] para τ≤t≤τ+ 65 indefinido para t > τ + 65 (2.13) donde Aes la constante de normalización relacionada con el pico de intensidad de corriente de la antena (en amperios), tes el tiempo en µs, y τes el conocido como envelope to cycle difference en µs. Idealmente, el envelope to cycle difference o ECD tiene un valor de τ= 0, por lo que se adoptará esta suposición a partir de ahora. Según [8], “El pulse trailing edge (aquella porción del pulso a continuación del pico del pulso, o a los 65 µs, cualquiera que se produzca primero) es controlada con el fin de mantener los requerimientos espectrales. Para diferentes lugares de transmisión, o con diferentes equipos de transmisión, el pulse trailing edge puede diferir significativamente en apariencia y características. A pesar de estas diferencias, para cada pulso y para todo t > 500 µs, i(t)satisface las tolerancias del pulse trailing edge, basadas en los picos de amplitud (A)” ([8, págs. 95-96], donde i(t)es lo que en este documento corresponde prácticamente con x0(t)). Estas tolerancias son: Categoría 1: i(t)≤0,0014A, (2.14) Categoría 2: i(t)≤0,016A. (2.15) Los requisitos espectrales se encuentran especificados en [13]: “la energía total fuera de la banda 90 −110 KHz debe de ser menor que un 1 % de la energía total radiada”. En el presente trabajo se ha optado por prolongar la forma de la envolvente definida para el intervalo τ≤t≤τ+ 65 para el caso t > τ + 65. Se utilizará esta simplificación a pesar de que no cumpla ni (2.14) ni (2.15) debido a que el establecimiento de estas tolerancias viene motivado únicamente por aspectos residuales de la señal que no intervienen de manera relevante en la extracción de información. El valor del parámetro A, se ha considerado darle el valor específico que haga que la señal completa se encuentre normalizada a 1, con el fin de facilitar posteriormente su interpretación. Para ello, se busca el valor de tque haga que E(t)sea máxima. Haciendo un análisis mediante el uso de las derivadas de la envolvente, se tiene que dE(t) dt = 2At exp −2t 65 1−t 65.(2.16) 12 Igualando a 0 y resolviendo para tse obtienen como posibles soluciones t= 0 yt= 65. Evaluando la segunda derivada para estos dos valores se obtiene d2E(t) dt2t=0 = 2A > 0 =⇒mínimo (2.17) d2E(t) dt2t=65 =−2Ae−2<0 =⇒máximo (2.18) con lo que es en el instante t= 65 µs cuando se alcanza el pico del pulso. Igualando a 1 el valor de la envolvente para ese instante, se obtiene el valor de: E(t= 65) = 652exp(−2·65/65)A=1=⇒A=e2/4225 ≃0,001749 (2.19) La envolvente E(t)modula a una portadora de manera que la señal modulada es x0(t) = E(t) sin(2πf0t)(2.20) donde f0es la frecuencia de la portadora, que en el sistema Loran-C es de 100 kHz. En lo sucesivo se tomará f0= 0,1MHz, con la única finalidad de que la escala temporal esté en µs. Este aspecto es sumamente sencillo de demostrar: f0[Hz]×t[s] = f0  [Hz]×1  106 [MHz] [ Hz]×t [s]× 106 1 [µs]  [s]=f0[MHz]×t[µs].(2.21) 0 50 100 150 200 250 300 350 400 450 t [ s] -1 -0.5 0 0.5 1 x0(t) Señal modulada Envolvente Figura 2.8: Forma del pulso unitario del sistema Loran-C. 2.2.2. Transmisión de los pulsos El patrón Loran-C más simple de transmisión es el mostrado en la Fig. 2.9, en el que varias estaciones slave, X, Y y Z, emiten pulsos individuales (como el mostrado en Fig. 2.8), mientras que la estación master, M, emite un doble pulso con el propósito de identificarla frente a las slaves. Si en vez de usar simple pulsing empleamos multiple pulsing, como en la Fig. 2.10, es posible incrementar la potencia media de 13 Figura 2.9: Transmisión Loran-C por single pulsing -desestimada- (Fuente: [14]). Figura 2.10: Transmisión Loran-C por multiple pulsing -real- (Fuente: [8]). los transmisores en un factor cuatro o mayor sin un correspondiente incremento de la potencia pico. Esta es la manera en la que se extiende al máximo la cobertura del sistema para una potencia pico dada [14]. Los pulsos desde el segundo al octavo de cada grupo se encuentran referenciados al primero de cada grupo, equiespaciados entre ellos 1000 µs más una cierta tolerancia (§4.3.5 en [8]). El noveno pulso correspondiente a la transmisión master está espaciado 2000 µs del octavo pulso del grupo. Este pulso se emplea primordialmente como ayuda a la identificación del grupo master y no como ayuda para la navegación [8]. El TD que ya se mencionó en §2.1.1 se esclarece gráficamente en Fig. 2.10, donde hay un TD por cada pareja master-slave. Nótese que las diferencias de tiempo se refieren siempre al primer pulso de cada grupo aunque se emplee multiple pulsing. Cada estación master del sistema Loran-C opera con un intervalo de repetición de grupo (más conocido como group repetition interval o GRI) que es múltiplo de 10 µs, desde 40000 µs hasta 99990 µs. Un GRI particular se identifica mediante el valor del GRI dividido entre 10, i.e. 7980 designaría un GRI de 79800 µs. El GRI no es más que el período entre dos recepciones consecutivas de pulsos correspondientes a un master. Los retrasos en las transmisiones de las estaciones secundarias con respecto a las master son seleccionadas de manera que se asegure el cumplimiento de los siguientes criterios para cualquier señal que pueda ser recibida de cualquier 14 chain particular [8]. La diferencia temporal mínima entre cualquier slave y el master es de 10900 µs. La diferencia mínima entre dos diferencias temporales es de 9900 µs. La diferencia temporal máxima es el valor del GRI menos 9900 µs. El espaciado mínimo entre los puntos correspondientes al último pulso del grupo de una estación y al primero del grupo consecutivo dentro del mismo chain es 2900 µs. El espaciado mínimo entre el noveno pulso del grupo master y el siguiente pulso slave (del mismo chain), sin embargo, debe de ser de tan sólo 1900 µs. Este es un resultado directo de la aplicación de los tres primeros criterios. Figura 2.11: Restricciones para la asignación del ED (Fuente: [8]). Cuando se introdujo en §2.1.1 el parámetro CD no se especificó ninguna limitación al valor que pudiera tomar, salvo ser mayor o igual que cero, ya que un retraso negativo no tendría significado físico. Con la necesidad del cumplimiento de estas restricciones en TD, CD se encuentra ahora acotado. Esto es debido a que, tal y como se mostró en (2.9), los valores extremos que puede tomar TD dependen de CD y BTT. Como este último únicamente depende de la posición espacial de las estaciones emisoras y de la velocidad de propagación de la señal (ambas supuestas constantes), CD será la única variable que permita ajustarse para satisfacer dichas restricciones. Recuérdese que existirá un valor de cada CD por cada estación slave del chain. Haciendo uso de lo anteriormente explicado podemos partir de (2.20) para continuar el desarrollo matemático hasta conseguir una expresión de la señal global. Los trenes de pulsos se pueden obtener mediante la convolución del pulso individual con un conjunto de impulsos [15]: xs1(t) = x0(t)∗Is(t) xm1(t) = x0(t)∗Im(t))para 0≤t < GRI (2.22) 15 donde ∗denota el producto de convolución2, e Is,m(t)son funciones que dependen de si los pulsos individuales x0(t)pertenecen a una emisión proveniente de una estación slave o de una estación master: Is(t) = 8 X m=1 δ(t−TD −1000(m−1)) ·pcsA(m)+ +δ(t−TD −GRI −1000(m−1)) ·pcsB(m)(2.23) Im(t) = 8 X m=1 δ(t−1000(m−1)) ·pcmA(m)+ +δ(t−GRI −1000(m−1)) ·pcmB(m) +δ(t−9000) ·pcmA(9) + δ(t−GRI −9000) ·pcmB(9)(2.24) donde δ(t)representa la función delta de Dirac, y las funciones pcXY son los parámetros de fase (phase-code parameter), que toman los valores mostrados en la Tabla 2.1 en función del carácter X (m -master- o s -slave-) e Y (A o B). Para su identificación, el primer grupo de pulsos en la secuencia se etiqueta como grupo A y el segundo grupo, un GRI después, es etiquetado como grupo B. Una secuencia de transmisión (phase-code interval o PCI) completa consta de un grupo A y un grupo B. A partir de ahí, la secuencia PCI se repite cada 2·GRI µs [8]. Phase Coding m 123456789 Master A +1 +1 -1 -1 +1 -1 +1 -1 +1 B +1 -1 -1 +1 +1 +1 +1 +1 -1 Slave A +1 +1 +1 +1 +1 -1 -1 +1 ∗ B +1 -1 +1 -1 +1 +1 -1 -1 ∗ Tabla 2.1: Valores que adoptan las funciones pcXY. El uso de la modulación de fase atiende a diferentes propósitos [14]: Aporta un método de identificación de las señales master yslave para los equipos automáticos (que sintonizan los chains de manera autónoma). Aceleran el proceso de búsqueda automática y simplifican los equipos receptores. La codificación de fase resuelve el problema de multihop skywave interference (Véase [14]). El método de codificación de fase proporciona comunicaciones entre estaciones Loran-C para propósitos administrativos. Asiste al receptor a operar cuando las señales son débiles en presencia de ruido [8]. 2Conviene recordar que el operador matemático de la convolución de dos funciones fygse define mediante f(t)∗g(t) := R∞ −∞ f(η)·g(t−η)dη. 16 Usando las propiedades de la convolución aplicadas a la delta de Dirac f∗(ag) = a(f∗g)∀a∈C(2.25) f(t)∗δ(t−t0) = f(t−t0)(2.26) en (2.22)-(2.23) y (2.22)-(2.24), se tiene xs1(t) = 8 X m=1 x0(t−TD −1000(m−1)) ·pcsA(m)+ +x0(t−TD −GRI −1000(m−1)) ·pcsB(m)(2.27) para el caso slave, y xm1(t) = 8 X m=1 x0(t−1000(m−1)) ·pcmA(m)+ +x0(t−GRI −1000(m−1)) ·pcmB(m) +x0(t−9000) ·pcmA(9) + x0(t−GRI −9000) ·pcmB(9)(2.28) para el caso master. Como la secuencia de transmisión se repite cada 2·GRI µs: xs,m(t) = xs1,m1(mod(t,2·GRI)),∀t≥0,(2.29) en el que mod(∗,∗) representa la operación módulo3. Finalmente la señal total que llega al receptor será la superposición de las señales enviadas por cada una de las estaciones del chain. x(t) = xm(t) + n X i=1 xs(t).(2.30) La señal completa dependerá del número de estaciones slaves (n) que dispone el chain en cuestión (n= 3 en este trabajo). El efecto del canal hace que a (2.30) se le agreguen señales espurias que la corrompan y degraden su calidad. 2.3. Propagación 2.3.1. Mecanismos de propagación El suelo, la troposfera y la ionosfera son responsables de que el modelo ideal de propagación en espacio libre no sea totalmente fidedigno en la mayoría de los casos reales. La orografía del suelo y sus características morfológicas, que condicionan sus propiedades eléctricas, afectan a la propagación de las ondas electromagnéticas. La banda de frecuencia de trabajo va a definir los mecanismos de propagación de ondas que hay que considerar a la hora de analizar un canal de radio. La señal Loran-C se encuentra en la banda de frecuencias LF en la que se produce la propagación principalmente mediante ondas de tierra (ground waves). Esta onda se propaga en la discontinuidad tierra – aire debido a las corrientes inducidas en la Tierra. El 3La operación módulo mod(a,n)devuelve el resto de la división euclídea de un número apor otro n. Esta operación matemática es útil a la hora de repetir una cierta secuencia numérica cada cierto período, en este caso n. 17 alcance que se obtiene varía con la frecuencia, la potencia transmitida y el tipo de suelo (tierra seca, húmeda, mar...) [16]. A su vez, para estas frecuencias también existe otros mecanismo de propagación denominado reflexión ionosférica (sky waves), en el que la presencia de la ionosfera hace que se reflejen las ondas radioeléctricas, provocando que éstas retornen a la Tierra. La ionosfera es la región superior de la atmósfera (aproximadamente entre los 60 a 400 km de altura) compuesta por una mezcla de gases ionizados (plasma) impenetrable para las ondas de radio de frecuencias por debajo de los 30 MHz, lo que provoca esta reflexión. Esta ionización, o presencia de electrones libres, se produce fundamentalmente por la absorción de radiaciones solares altamente energéticas, como las radiaciones ultravioleta, rayos X y radiación gamma. Esto hace que la densidad de electrones varíe según la hora del día, la estación del año y los ciclos de manchas solares. La ionosfera se encuentra dispuesta en las siguientes capas [16]: La primera capa es la capa D, producida entre los 60 y los 90 km. Decrece rápidamente tras la puesta de sol, por lo que sólo aparece durante el día. Es sumamente absorbente para frecuencias por debajo de unos 10 MHz, ocasionando que las señales lleguen a la Tierra muy atenuadas. La capa E está localizada entre los 90 y 130 km. Aparece fundamentalmente de día, y muy tenuemente por la noche. La capa más importante, donde se produce principalmente la reflexión ionosférica, es la capa F. Se extiende desde los 150 hasta los 400 km, y de día se desdobla de dos capas: F1 y F2. La desaparición de la capa D y la permanencia de la capa F durante la noche origina que la cobertura del sistema Loran-C sea mayor en este período que durante el día (véase Fig. 2.13), ya que la desaparición de la capa absorbente permite a las ondas llegar hasta la capa F, reflejarse y alcanzar lugares hasta tres veces más distantes de la Tierra que las ground waves. El área de cobertura se establece como aquella zona en la que la recepción de la señal tiene un nivel que permita determinar la posición con un grado de exactitud y SNR específico [8]. La operación normal de Loran-C asume la recepción de ground waves para obtener una localización con una alta precisión, ya que son señales más estables [14]. Las sky waves siempre llegan después que las ground waves, aunque la diferencia entre ambas se hace cada vez más pequeña a medida que aumenta la distancia al transmisor, no siendo nunca menor que 30 µs en ningún punto del área de cobertura. Esta condición asegura por tanto que durante los primeros 30 µs de un pulso Loran-C no se producirá una contaminación por sky wave. Para distancias superiores a las 1000 NM4, las ground waves es probable que sean inutilizables debido a que sufren una atenuación mayor que las sky waves. Por consiguiente, las sky waves se emplean para alcances mayores que este valor, siendo consciente de que estas ofrecen una menor precisión y que es necesario aplicar factores correctivos para compensar la diferencia en la trayectoria seguida por la onda en comparación con la ground wave [8]. Por otro lado, se ha constatado que la fiabilidad del sistema Loran-C se ve sólo ligeramente afectada por las condiciones atmosféricas adversas, a diferencia de otros sistemas de radionavegación. 4Nautical miles 18 Figura 2.12: Ground waves ysky waves (Fuente: [11]). 2.3.2. Velocidad de propagación La velocidad de propagación de una señal es v=c0 √r ,(2.31) donde c0es la velocidad de la luz en el vacío y res la permitividad relativa (o constante dieléctrica) del medio a través del cual se propaga dicha señal. La permitividad relativa depende esencialmente del material, y en menor grado de otros factores como la frecuencia, la temperatura, etc. Estos factores provocan que cualquier medida de tiempo pueda estar ligeramente equivocada. Con el objetivo de que el sistema opere dentro de los límites de exactitud que se requiere, se deberán aplicar unos factores de corrección. Hay tres factores de corrección que pueden ser aplicados y que se denominan factores de fase [8]. Factor de fase primario (PF): relacionado con el hecho de que la velocidad de propagación atmosférica de la señal es ligeramente más lenta que en el vacío. Esta diferencia es debida a que el índice de refracción de la atmósfera es ligeramente superior a la unidad (r≥1). Factor de fase secundario (SF): relacionado con el hecho de que la velocidad de propagación de la señal es más lenta cuando viaje por el agua marina ya que su conductividad es mayor que comparada con la tierra. Este factor permite corregir el tiempo extra necesario para la propagación de la señal por un camino completo de agua marina comparado con uno sobre tierra. Factor de fase secundario adicional (ASF): debido a que los transmisores Loran- C se encuentran situados en tierra, la señal propagada viajará sobre tierra y agua. El ASF puede ser calculado tratando el camino de propagación como segmentos separados, cada uno de los cuales posee un valor de conductividad dependiente de si el segmento es sobre tierra o sobre agua. Además, factores 19 PE: elemento que computa la estimación de la posición. En caso de presentarse ambigüedades también es el encargado de solventarlas. 3.2. Detección de señal 3.2.1. Filtros (BPF y LPF) El ruido y las interferencias son, en prácticamente la mayoría de las aplicaciones, indeseables. Por esta razón, en multitud de sistemas de comunicaciones se emplean técnicas para su supresión o atenuación para así obtener una señal de mejor calidad comparada con la original. La herramienta aquí empleada para cumplir este propósito será el uso de filtros FIR (finite impulse response), designados de esta manera por disponer de una respuesta impulsional de duración finita. Para un filtro FIR causal en tiempo discreto de orden N, cada valor de la secuencia de salida es la suma ponderada de los valores más recientes de la entrada al filtro [18]: y[n] = N X k=0 bk·x[n−k](3.1) donde x[n]es la señal de entrada, y[n]es la señal de salida, Nes el orden del filtro (un filtro de orden N-ésimo tiene N+1 coeficientes bk), y bkes el valor de la respuesta impulsional en el instante k-ésimo para 0≤n≤Nen un filtro de orden N-ésimo. Figura 3.2: Forma directa de un filtro FIR discreto de orden N. Cada unidad de retraso es el operador z−1en notación de la transformada Z (Fuente: [18]). En concreto, se emplearán dos filtros FIR de orden 32: el primero será un filtro con comportamiento BPF (band pass filter) con un ancho de banda BW centrado en f= 100 kHz para reducir el nivel de ruido inmediatamente después de recibir la señal Rx, y en segundo lugar un filtro con comportamiento LPF (low pass filter) con frecuencia de corte f=BW/2situado a la salida del DdE con el propósito de disminuir el rizado de la envolvente detectada. Las respuestas en frecuencia de ambos filtros utilizados se muestran en la Fig. 3.3. El correcto comportamiento de los filtros ante las señales que atañen a este receptor Loran-C se comprueba con las respuestas en frecuencia de la Fig. 3.4. Cada una de las gráficas de esta figura corresponde a un filtrado determinado: la superior al de paso banda y el inferior al de paso bajo. En ambas situaciones el proceso de filtrado pasa una señal input (color azul) a una output (color rojo) con la intención de 26 0 50 100 150 200 250 300 350 400 450 500 Frecuencia [kHz] -100 -50 0 Módulo (dB) Filtro BPF 0 50 100 150 200 250 300 350 400 450 500 Frecuencia [kHz] -100 -50 0 Módulo (dB) Filtro LPF Figura 3.3: Módulo de la respuesta en frecuencia de los filtros FIR BPF y LPF. asemejarla lo máximo posible a la ideal (color verde), que corresponde a la respuesta de las señales en el transmisor. Todo filtrado de una señal introduce un desfase. Esto significa que la señal de salida del filtro se encontrará desplazada en tiempo con respecto a la de entrada. Los filtros de respuesta impulsional finita a menudos provocan un retraso de igual cantidad para todas las componentes de frecuencia, lo que los hace fáciles de compensar en términos temporales [20]. Se puede demostrar que el retraso de estos filtros es de τ=N/2(cada filtro). 3.2.2. Detector de envolvente El detector de envolvente (DdE) digital se corresponde con el circuito mostrado en la Fig. 3.5. La diferencia fundamental entre este DdE con su equivalente analógico es la sustitución del diodo por un interruptor dinámico controlado, de tal manera que el circuito se encontrará en modo ON (dinámica rápida o de carga) cuando dicho interruptor se encuentre conectado y en modo OFF (dinámica lenta o de descarga) cuando esté desconectado. Para ambos casos, la relación entre el voltaje de entrada x(t)y el de salida y(t)se detallan a continuación. 27 0 200 400 600 800 1000 Frequency (kHz) -180 -160 -140 -120 -100 -80 -60 -40 PSD (dB/Hz) PSD de señales en etapa BPF Tx Rx Rx BPF 0 200 400 600 800 1000 Frequency (kHz) -160 -140 -120 -100 -80 -60 -40 -20 PSD (dB/Hz) PSD de señales en etapa LPF env Tx env Rx env LPF Figura 3.4: PSD (Power Spectral Density) de las señales entrante y saliente a los filtros BPF y LPF, así como el de su original transmitida. Modo ON Cuando el interruptor se encuentra conectado, la ec. que gobierna la dinámica del circuito es x(t) = R2y(t) R1 +Cdy(t) dt +y(t).(3.2) 28 Figura 3.5: Circuito DdE digital. Discretizando el dominio temporal5para la ec. (3.2) y haciendo uso de las relaciones K1=R1C Ts =τ1 Ts⇒dinámica lenta,(3.3) K2=R2C Ts =τ2 Ts⇒dinámica rápida,(3.4) se tiene que x[n] = K2 K1 y[n] + K2(y[n]−y[n−1]) + y[n],(3.5) donde despejando la salida y[n] = K1 K1+K2+K1K2 x[n] + K1K2 K1+K2+K1K2 y[n−1].(3.6) Asumiendo que K1K2, podemos simplificar la ec. (3.6) para obtener y[n] = 1 1 + K2 x[n] + K2 1 + K2 y[n−1].(3.7) En (3.7) puede verse que el primer término corresponde a la actualización y el segundo a la memoria. Modo OFF Cuando el interruptor se encuentra desconectado sólo entran en juego el condensador Cy la resistencia R1. La descarga del condensador proporciona una salida no dependiente de la entrada x(t)del circuito. y(t) = −R1Cdy(t) dt .(3.8) Procediendo de manera análoga al apartado anterior se llega a y[n] = −K1(y[n]−y[n−1]) = K1 1 + K1 y[n−1].(3.9) Debido a que en esta ocasión el circuito se encuentra aislado de la entrada, sólo se tiene término de memoria. 5El proceso de discretización cumple que x[n] = x(t)|t=nTs=x(nTs)∀n∈N0. Para la discretización de las derivadas se ha empleado la diferencia dividida regresiva de orden uno. 29 Reuniendo ambas coyunturas se tiene que a la salida del DdE será y[n] =        1 1 + K2 x[n] + K2 1 + K2 y[n−1] =⇒modo ON K1 1 + K1 y[n−1] =⇒modo OFF (3.10) En esta ocasión, ya que el control del interruptor es digital, hará falta establecer una condición que haga que el circuito permute entre los estados ON y OFF. Esta condición no es más que (modo ON ⇐⇒ y[n−1] ≤x[n], modo OFF ⇐⇒ y[n−1] > x[n],(3.11) y es simple de comprender: cuando el valor de la señal modulada x[n]alcanza el valor de la envolvente y[n], el circuito se establece en modo ON y la envolvente la ‘persigue’ durante su ascensión hasta el momento en el que x[n]comienza a decrecer de manera más acusada de la que lo hace la envolvente, instante en el que se permuta al estado OFF. Los parámetros de ajuste del DdE son τ1yτ2, por lo que el desempeño de este dependerá de los valores particulares que se le asignen. A modo ilustrativo, en la Fig. 3.6 se presenta el comportamiento del DdE para diferentes valores de estos parámetros. Ambas gráficas de la fila superior corresponden a cambios en τ1por encima y por debajo del valor que finalmente se ha ajustado, mientras que las inferiores corresponden a cambios en τ2, igualmente por encima y por debajo del valor ajustado. Con ellas se puede comprobar que: A mayor valor de τ1, al estar relacionado con la dinámica lenta, la descarga será más lenta y a la inversa. A menor valor de τ2, al estar relacionado con la dinámica rápida, la carga será más rápida y a la inversa. Finalmente se han ajustado estos parámetros a los valores τ1= 50 µs y τ2= 0,5µs, que proporcionan la respuesta que se muestra en la Fig. 3.7. Es importante destacar que el empleo de un DdE implica que la información del phase coding se abandone ya que se pasa a trabajar con la envolvente en vez de con la señal oscilatoria modulada, que es donde va el PC. El DdE es, en definitiva, un caso particular de filtro no lineal, lo que dificulta el conocimiento analítico del retraso que introduce. No obstante, se ha podido estimar que para los valores τ1yτ2ajustados se tiene una demora en el entorno a las 9 muestras. Si aumentara (disminuyera) el valor de alguno (o ambos) de los parámetros, el valor del retraso aumentaría (disminuiría). 3.2.3. Correlación cruzada En la disciplina del procesamiento de señales la correlación cruzada o crosscorrelation (CC) es una medida de la similitud entre dos señales, frecuentemente usada para encontrar características relevantes en una señal desconocida por medio 30 0 100 200 300 400 -0.5 0 0.5 x(t) 1 = 100 s 2 = 0.5 s Señal Envolvente 0 100 200 300 400 -0.5 0 0.5 1 = 25 s 2 = 0.5 s Señal Envolvente 0 100 200 300 400 t [ s] -0.5 0 0.5 x(t) 1 = 50 s 2 = 5 s Señal Envolvente 0 100 200 300 400 t [ s] -0.5 0 0.5 1 = 50 s 2 = 0.025 s Señal Envolvente Figura 3.6: Diferentes comportamientos del DdE ante distintos valores de τ1yτ2. de la comparación con otra que sí se conoce (patrón o referencia). Para el caso particular de funciones discretas reales, la correlación cruzada entre dos señales fy gse define tal que: f[n]? g[n] := X k f[k]·g[n+k],(3.12) donde la sumatoria se realiza sobre valores enteros de kapropiados.6La correlación cruzada es útil para determinar el retraso temporal entre dos señales. Una vez calculada la CC de ambas, el máximo de la función de CC indicará el instante temporal para el cual las señales se encuentran lo más alineadas posible [18], pudiendo así conocer el retraso entre ellas. Esta herramienta se aplicará de la siguiente manera: la señal env_LPF (envolvente filtrada, Fig. 3.1) se correlará con una señal patrón (denominada REF_PULSE) que será la envolvente ideal de un pulso unitario Loran-C (en Fig. 2.8, línea sólida roja en intervalo t∈[50,400]). El resultado de esta operación puede contemplarse en la gráfica inferior de la Fig. 3.14, en la que puede comprobarse que los máximos de CC coinciden con los instantes en los que se inician cada uno de los pulsos individuales de la señal (TOA, time of arrival) (gráfica superior de la misma Fig.). La robustez de este instrumento reside en el calificativo de ruido blanco, que signi- 6La correlación cruzada se relaciona con la convolución mediante f(t)? g(t) = f∗(−t)∗g(t). 31 0 50 100 150 200 250 300 350 400 450 t [ s] -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 x(t) 1 = 50 s 2 = 0.5 s Señal Envolvente Figura 3.7: Comportamiento del DdE para los parámetros τ1yτ2ajustados. fica que sus valores no se encuentran correlacionados entre sí. Esto permite que al correlar nuestra señal con ruido con un patrón limpio, los instantes de alineación puedan ser fácilmente identificables. El offset que aparece en la función CC es debida a esta característica del ruido AWGN. 3.3. Control del receptor Las actuaciones concernientes a esta etapa del procesado vienen gobernadas mediante una señal de control biestado denominada DetectMode, la cual indicará si el receptor se encuentra en modo detección/búsqueda (MD) o en modo escucha (ME). Las conmutaciones entre ambos estados son realizadas de manera autónoma por el propio receptor, sin necesidad de manipulación externa, favoreciendo la automatización del sistema. Con este método se pretende que, durante aquellos intervalos temporales existentes entre dos conjuntos consecutivos de pulsos cualesquiera7(intervalo de ‘silencio’), el receptor se encuentre midiendo el ruido adherido a la señal por efecto de su paso por el canal (ME); mientras que cuando se detecte que se está recibiendo uno de estos conjuntos, el receptor se encargará de determinar los instantes en los que se inician cada uno de los pulsos que lo conforman, así como identificando cada conjunto como master oslave y su posterior posprocesado hasta conseguir la estimación de la posición (MD). Así mismo, desde el punto de vista de la implementación real en hardware del receptor, donde las restricciones de diseño suelen venir dadas principalmente por la capacidad computacional y de memoria de 7Independientemente de la estación emisora del que provengan. 32 los componentes, la existencia de un modo de operación ‘relajado’ reduce la carga de trabajo del microprocesador, lo que genera una ventaja significativa. No se pretende en esta Memoria explicar de manera pormenorizada el código desarrollado en matlab para el receptor Loran-C, sino una descripción sucinta y concisa del funcionamiento del mismo. 6.5 7 7.5 8 8.5 9 9.5 10 10.5 11 11.5 t [ s] 104 0 20 40 60 80 100 120 140 160 180 CC(t) CC DetectMode CheckMax RAW_DTH CC_DTH noise_floor noise_floor+noise_std Figura 3.8: Umbrales, medidas del ruido e instantes CheckMax. En la Fig. 3.8 se muestra un ejemplo ilustrativo de esta señal de control (línea sólida roja), en la que un valor ‘alto’ corresponde al MD y un valor ‘bajo’ al ME. La referencia del valor bajo de esta señal se ha tomado caprichosamente a una altura arbitraria con la única finalidad de facilitar la visualización simultánea de la señal de control con la CC. 3.3.1. Modo escucha 3.3.1.1. Activación del modo escucha La transición de MD a ME acontecerá cuando, una vez se produzca el inicio del MD, haya transcurrido una ventana temporal, denominada RAW_DETECTION_WINDOW, lo suficientemente grande como para abarcar la duración del conjunto de pulsos más largo existente en el sistema Loran-C, asegurando así que durante la recepción de pulsos el sistema se encuentra percibiéndolos. Se recuerda que el sistema Loran-C presenta dos configuraciones de conjuntos de pulsos, los relativos a una estación slave y los correspondientes a una estación master. Los primeros se componen de 8 pulsos equiespaciados una cantidad de 1000 µs (denominada W_1), mientras que los segundos tienen la misma estructura que los ya mencionados, con la salvedad de que se le anexa un noveno pulso espaciado con su inmediatamente precedente el doble de W_1. Con ello, es fácil comprobar que el conjunto de pulsos más largo es el 33 relativo al master, con una duración máxima8de 8×W_1 + 2 ×W_1 = 10000 µs, que será el valor adoptado por RAW_DETECTION_WINDOW. 3.3.1.2. Actuaciones durante el modo escucha Como se ha admitido que las interferencias que afectan a la señal se pueden modelar mediante un ruido blanco gaussiano aditivo (AWGN), la señal recibida entre conjuntos de pulsos se trata exclusivamente del ruido del canal. Por tanto, mientras el modo de escucha se encuentra activado, el receptor únicamente se encarga de realizar una medidas sobre el ruido; en particular de su media y su desviación estándar. Estos valores, almacenados en las variables noise_floor ynoise_std respectivamente (matemáticamente µyσ), se reservarán para establecer unos umbrales que se explicarán a continuación. Ambas medidas estadísticas requieren disponer de un conjunto de valores del ruido anteriores al instante actual. Para ello, se toma una ventana de duración INACTIVE_TH que va guardando los últimos valores anteriormente medidos del ruido. La elección de la longitud de esta ventana debe ser un compromiso entre: un vector lo suficientemente largo para que las medidas estadísticas estén convenientemente enriquecidas, pero no de manera excesiva como para no permitir ni una sola estimación de estos valores. Los umbrales que se calculan a partir de las medidas del ruido son: RAW_DTH = (µ+σ)∗ρ+ 1,(3.13) CC_DTH =RAW_DTH ∗η, (3.14) donde RAW_DTH se define como el valor que debe superar la señal CC para pasar de ser considerada ruido a ser considerada como ‘señal útil’, y donde CC_DTH se establece como el límite por debajo del cual el valor que se encuentre almacenado como máximo local actual no será consolidado como máximo9(ambos umbrales se muestran en un ejemplo de señal en la Fig. 3.8). En estas ecuaciones µrepresenta la media, σla desviación estándar, y ηun coeficiente que será explicado en detalle a continuación. No se debe caer en la tentación de presuponer que los valores de estos umbrales pueden ser elegidos de manera caprichosa; es más, son una de las piezas clave del buen o mal desempeño del receptor. Es por eso que se hace se torna imperativo el detenerse un instante para justificar las ecs. (3.13) y (3.14). En primer lugar se calcula la media estadística del ruido en una ventana muestreal de longitud INACTIVE_TH, que nos proporcionará una estimación del nivel de suelo en el que se encuentra centrado (noise_floor en Fig. 3.8), a la vez que se obtiene también la desviación estándar de la misma ventana. Este valor dará una medida del grado de dispersión de los datos con respecto al valor promedio. Gracias a que el ruido se comporta como una variable aleatoria con distribución normal, el conocimiento de estas cantidades permite predecir un rango de valores dentro del cual la variable aleatoria tomará valor con una probabilidad determinada, tal y como se muestra en la Fig. 3.9 para el caso teórico gaussiano. Se puede observar en ella que a medida que el rango de valores [µ−%σ, µ +%σ], centrado en µ, aumenta en longitud (aumentando %), la probabilidad de contener en él un valor tomado por la variable aleatoria aumenta. En particular, y debido a que 8Se considera que los pulsos Loran-C son prácticamente cero para valores temporales de más de 1000 µs tras el inicio de este. 9La comprensión de estos umbrales se verá reforzada en los apartados posteriores. 34 Figura 3.9: Porcentaje de los valores aleatorios que caen dentro de la banda centrada en la media de una distribución normal y con un ancho de dos, cuatro y seis desviaciones estándar (Fuente: [18]). se está interesado en establecer un umbral por encima de la media, sólo nos atañerá el límite superior del intervalo indicado. En esta ocasión, tras una proceso de ajuste empírico tras diversas simulaciones, se ha optado por tomar un umbral superior de la forma ρ(µ+σ), donde ρ= 1,4para asegurar que se encuentra lo suficientemente apartado verticalmente del ruido como para que los valores de la señal que alcancen este límite puedan ser tomados como ‘señal útil’ con una certeza de ≈83,85 %. Esta probabilidad se calcula haciendo uso de la teoría de probabilidad, que dicta que P(ρ(µ−σ)≤X≤ρ(µ+σ)) = Φµ,σ2[ρ(µ−σ)] −Φµ,σ2[ρ(µ+σ)] (3.15) donde Xes la variable aleatoria ruido y Φµ,σ2(x)es la función de distribución de la distribución normal definida como sigue: Φµ,σ2(x) = 1 σ√2πZx −∞ exp −(u−µ)2 2σ2du (3.16) Tomando para este cálculo aproximado la función de distribución de la normal estándar, donde µ= 0 yσ= 1, y particularizando en (3.15), la probabilidad que se obtiene es de ≈83,85 %. Si se quisiera un valor de certeza mayor, no habría más que aumentar el valor de ρ. Se puede demostrar fácilmente que los coeficientes ρy %son intercambiables, siendo su relación %=ρ(µ+σ)−µ σ.(3.17) En la ec. (3.13) aparece un segundo sumando de valor unidad. Este offset se implanta porque en caso de que no existiera ruido en canal, µ=σ= 0, lo que supondría que RAW_DTH=0. Esta circunstancia provoca que ∀CC(t)>0cuando el ME está activado, aún a pesar de que debiera teóricamente debiera tomar 010, activaría el MD innecesariamente. 10Coyuntura habitual cuando se trabaja con un sistema numérico de precisión finita. 35 2 4 6 8 10 12 14 16 t [ s] 104 -1.5 -1 -0.5 0 0.5 1 1.5 x(t) Rx Rx BPF Tx env LPF 2 4 6 8 10 12 14 16 t [ s] 104 0 50 100 150 200 CC(t) CC Slave Master DetectMode RAW_DTH CC_DTH Figura 3.14: Señales producidas de cada etapa del receptor. 42 3.4.1. Fuentes de error 3.4.1.1. Dependencia de TD con el tiempo Todo sistema de posicionamiento basado en la diferencia de dos (o más) medidas de tiempos posee un error inherente al propio método de posicionamiento. Este hecho se ilustra en la figura a continuación: Figura 3.15: Error de estimación de la posición debida al movimiento del vehículo. Si se supone que el sistema se halla embarcado en un vehículo que, durante el intervalo de tiempo ∆tentre la primera medida de tiempo (tX) y la segunda (tY> tX), se encuentra en movimiento, es evidente que en ambos instantes el vehículo se encontrará en posiciones del plano P(tX)yP(tY)diferentes por el mero hecho de estar desplazándose. Esta situación hace que la pareja de hipérbolas que se toman para intersectarlas (ΓX(tX)∩ΓY(tY)) y obtener así una estimación de la posición ˜ P(que pudiera darse en un cierto instante posterior11,tY+τ, a la segunda medida temporal), correspondan a posiciones diferentes. Esto hace que, como bien se ilustra en la Fig. 3.15, la posición estimada ˜ P(tY+τ)y la posición real P(tY+τ)difieran en una cantidad (segmento de color cian) que será significativa o no dependiendo de la velocidad del vehículo, ya que a mayor velocidad, mayor distancia recorrida durante ∆ty por tanto mayor error en la estimación. Nótese que en el caso particular del reposo, ambas hipérbolas sí corresponderían a la misma posición, ya que durante el intervalo ∆tel vehículo permanece fijo si no existe movimiento, por lo que el error cometido entre la posición estimada y la real desaparecería. La conclusión de la coyuntura planteada es que, siempre que el sistema de posicionamiento Loran-C se encuentre en un vehículo en movimiento, ello implica que: Existirá un error ineludible entre la posición estimada y la real que será significativo o no dependiendo de la velocidad del vehículo. 11En el código de matlab, este instante se denomina PRT (position refresh time). 43 Ya que la posición variará con el tiempo, el vector TD pasará a ser dependiente de este también (TD =TD(t)), como ya se había indicado. 3.4.1.2. Resolución del receptor Como todo sistema digital, el receptor operará a una cierta tasa o período de muestreo Tsdurante el proceso de discretización de la señal analógica recibida (sampling), siendo este el tiempo entre intervalos de toma de muestras. En lo sucesivo, se denominará resolución del receptor al menor incremento con el que el receptor es capaz de obtener una medida, que coincide con Ts. Esta discretización del tiempo limita el desempeño del receptor en la búsqueda de los instantes de inicio de los trenes de pulsos (tMytSi, con i=X,Y,Z...), ya que ∀t≥0tal que mod(t,Ts)6= 0 (es decir, para todo instante de tiempo distinto de los múltiplos de Ts), el receptor sólo será capaz de proporcionar una estimación del t real de manera que ˜ t=bt/TscTs∈0,Ts,2Ts,...nTs∀n∈N0, y donde b∗cexpresa la función parte entera. Esto ocasiona un error de resolución temporal tal que δt < Ts, que se propagará al cálculo del TD y que se denotará como δTD. 9 9.5 10 10.5 11 11.5 12 12.5 t [ s] 0 0.5 1 1.5 2 2.5 3 3.5 4 x(t) 10-3 Ts = 3 s Ts = 1 s Ts = 0.01 s Inicio pulso Figura 3.16: Señal formada en el receptor para diferentes valores de Ts. Puede verse que cuando el instante de inicio del pulso no es un múltiplo entero del Tsseleccionado, el valor estimado del inicio diferirá del real en una cantidad δt tanto mayor cuanto mayor sea el valor de Ts. 3.4.1.3. Gradiente hiperbólico Lo que se pretende a continuación es caracterizar el error de posicionamiento que se comete si el receptor, en vez de posicionarnos sobre una LOP correspondiente a un cierto TD, lo hiciera sobre una que correspondiera a TD +δTD. Desde un punto de vista pragmático, al operador del receptor lo que le interesará conocer será el área de incertidumbre de la localización estimada. Este área es entendida como aquel sector situado en los alrededores de la posición estimada donde 44 cambios en el TD menores a la resolución del receptor δTD no serán detectados por este. Para ilustrar esta idea véase la Fig. 3.17. Figura 3.17: Áreas de incertidumbre en función del gradiente hiperbólico. En esta figura se han representado una serie de LOPs arbitrarias (líneas sólidas) correspondientes a dos parejas master-slave con master común, las LOPs que resultan de sumar o restar al TD de la principal la cantidad δTD (líneas punteadas) y las áreas de incertidumbre (área sombreadas). Las líneas de posición LOP2yLOP4 pertenecerán a una de las parejas, mientras que LOP1yLOP3corresponderán a la otra. Según las LOPs particulares que estime el receptor y su posterior intersección, hallaremos la posición del vehículo. P1=LOP1∩LOP2,(3.22) P2=LOP1∩LOP4,(3.23) P3=LOP3∩LOP4.(3.24) Puesto que el receptor sólo podrá dar estimaciones discretas de TD (... TD − δTD, TD, TD +δTD ...), sus LOPs correspondientes tendrán unos pasos fijos entre sus adyacentes. Por lo tanto, hay que destacar dos cuestiones: Un familia discreta de LOPs que se intersecan con otra familia diferente da lugar a un mallado de puntos en el plano que hará que, toda aquella posición real del vehículo no coincidente con un punto de la malla, se estime por uno que si pertenece a ella. Por otro lado, los errores en las estimaciones del TD por parte del receptor, que siempre que se produzcan serán como mínimo de ±δTD, harán que en 45 vez de darse la posición ideal Pise dé como estimación alguno de los vértices del área sombreada. En caso de que el error cometido sea mayor, debido por ejemplo a una alta tasa de ruido, es evidente que vértices distarán una mayor cantidad de la posición ideal, incrementando el error de posicionamiento. Nótese que fijada la resolución del receptor, las áreas de incertidumbre se amplían a medida que el punto de intersección se aleja de la baseline, por lo que el gradiente entre las ramas hiperbólicas aumenta. En este contexto se define el gradiente hiperbólico como la distancia perpendicular a la tangente de un LOP en un punto que hay que recorrer para alcanzar el LOP consecutivo. Figura 3.18: Familia de LOPs de una pareja master-slave donde se ilustran los excesivos gradientes que se tienen cerca de la baseline extension (Fuente: [8]). Teniendo estas consideraciones presentes, ahora se está en condiciones de comprender el porqué de que, generalmente, los chains tengan un mayor número de estaciones que las estrictamente necesarias para el posicionamiento12. Para una localización particular dentro de la cobertura de más de un chain Loran-C, el operador 12Se recuerda que el mínimo son 2 estaciones slaves por cada master. 46 (de forma manual) o el propio receptor (de forma automática) debería seleccionar el mejor chain disponible, entendiéndose este como el que se reciba con mayor potencia13. Una vez se ha seleccionado el chain más conveniente, es necesario seleccionar la triada master más dos slaves que ofrezcan el mejor posicionamiento. Es aquí donde entra en juego la redundancia de estaciones slaves. Para una misma ubicación pero para distintas triadas, el área de incertidumbre puede diferir significativa entre ellas. Idealmente, dos LOPs que se corten en ángulo recto deberán ser siempre las que se utilicen, ya que son las que darán mejor precisión. Dado que esto no es siempre posible de conseguir, las LOPs que se deben elegir serán las que intersequen lo más cerca posible al ángulo recto. Este requisito es equivalente a buscar la triada con menor área de incertidumbre, o lo que es lo mismo, menor gradiente entre LOPs. En la regiones cercanas al baseline extension de una pareja master-slave, los gradientes se hacen muy grandes, potenciando que se produzcan errores muy grandes. Cuando esto ocurre, el sistema Loran-C provee de un patrón repetitivo de parpadeo on-off (aproximadamente 0,25 s ‘on’ y 3,75 s ‘off’) de los dos primeros pulsos de la señal slave, el cual advierte que esa estación es inservible en la configuración actual (para más detalle, véase §4.3.6 en [8]). 3.4.2. Fuentes adicionales Man-made noises [23]: •Cross-rate interference (CRI): provocado por aquellas señales originadas en un chain diferente que opera en las cercanías del que se ha sintonizado. El problema proviene por la corrupción que puede provocar el hecho de que todas las estaciones empleen una misma frecuencia de transmisión. •Carrier-wave interference (CWI): el espectro de Loran-C va desde los 90 hasta los 110 kHz. Otras señales que se encuentren en el espectro circundante de los 50-150 kHz son potencialmente perjudiciales para la operación de Loran-C. Estas señales se las denominan CWI. Cambios en la velocidad de propagación de la señal de radio por cambios en el tipo de medio o condiciones ambientales en los mismos. Inestabilidades en el control del valor de CD. 13Es posible aplicar otro criterio de selección de chain si así se considera oportuno. 47 Capítulo 4 Resultados de las simulaciones 4.1. Metodología TRANSMISOR CANAL Definición trayectoria Cálculo TD Conformado señal Tx Comparación Análisis de prestaciones AWGN Tx RECEPTOR Figura 4.1: Esquema metodológico propuesto para el estudio del sistema Loran-C. El interior del bloque receptor es el exhibido en la Fig. 3.1. LaFig. 4.1 ofrece una visión general del planteamiento llevado a cabo para estudiar el sistema. Al igual que cualquier otro sistema de telecomunicaciones, este presenta tres bloques principales: transmisor, canal y receptor. La idea que se sugiere es la siguiente: en primer lugar predefinir una serie de puntos coordenados (trayectoria) que se supone serán recorridos por un vehículo, en este caso una aeronave. El conocimiento de esta trayectoria permite calcular la forma que tendrá la señal de radio (efectuada en el transmisor) a su llegada al receptor. Sin embargo, la existencia del canal entre Tx y Rx introducirá un ruido que perturbará la calidad de la señal. Estas interferencias harán que la información de posición que se extraiga de esta señal (efectuado en el receptor) difiera con respecto a la original. La cuestión es comparar la trayectoria definida inicialmente con la estimada por el receptor y analizar, en función de su similitud, el desempeño del sistema. 49 4.2. Simulaciones realizadas En este capítulo se muestran y comentan los resultados obtenidos de las diferentes simulaciones llevadas a cabo. Con ellos se trata de verificar el correcto funcionamien- to del código desarrollado para modelar el sistema Loran-C. Las simulaciones se clasificarán en escenarios diferentes: posicionamiento con el vehículo inmóvil, posicionamiento durante trayectoria recta y posicionamiento durante viraje. Se analizan estas tres configuraciones por ser las más habituales durante el vuelo de una aeronave. Para todos los casos la dinámica de estudio será la siguiente: para cada uno de los valores de SNR (Signal to Noise) en un rango que irá desde un valor alto (buena calidad de señal) hasta uno bajo (baja calidad de señal) se realizará una simulación lo suficientemente duradera como para poder obtener una buena estimación del error producido en el proceso de posicionamiento. En cada caso el error se caracterizará mediante el parámetro que se ha considerado más apropiado. Con la información suministrada por las simulaciones será posible confeccionar tanto una gráfica comparativa error-SNR como un mapa en coordenadas cartesianas planas donde se presenten y contrasten las posiciones estimadas frente a las reales. En las Tablas 4.1 y 4.2 se recopilan todos los parámetros que se han fijado para la realización de las diferentes simulaciones, posibilitando de esta manera una eventual reproducción de los resultados aquí obtenidos. Ocasionalmente algún parámetro muestra más de un valor, ya que en ciertas situaciones se han variado en un rango para comprobar distintos aspectos del comportamiento del receptor, como son Ts(y Fsconsecuentemente), la velocidad de la aeronave, el SNR y el tipo de trayectoria seguida. Cuando aparece el símbolo ∗en una casilla quiere decir que ese campo no procede rellenarlo. Todas las filas menos las pertenecientes a la categoría Loran_Rx son parámetros que se establecen en la función genConfig, mientras que el resto se establecen dentro del propio receptor. De los valores adoptados conviene reseñar que las frecuencias de muestreo cumplen, como no podía ser de otra forma, el teorema de Nyquist: fs>2·fmax (4.1) donde idealmente fmax = 100 kHz, la frecuencia de modulación. 4.2.1. Posicionamiento estático En caso de que se tenga la certeza de que el vehículo del que queremos determinar su ubicación se encuentre inmóvil, las LOPs sobre las que el sistema Loran-C lo posicionará serán las mismas (salvo errores de estimación de las medidas) para todo instante de tiempo; es decir, los valores de TD que identifican dichas LOPs no dependerán del tiempo. Antes de presentar la magnitud que se va a emplear para cuantificar el error, hay que clarificar una idea fundamental. Aunque en el lenguaje cotidiano los términos exactitud y precisión se usen indistintamente para referirse al mismo concepto, en las ciencias físicas y estadísticas es importante hacer una distinción entre ambos: la exactitud (accuracy) es la cuantificación de cómo de cerca se encuentran las medidas al valor “real” (valor aceptado como correcto), mientras que la precisión (precision) es la cuantificación de cómo de cerca se encuentran las diferentes medidas unas de otras [24]. Figuras como 4.2 son comúnmente utilizadas para esclarecer las diferencias 50 Parámetro Descripción Valor Unidad Ctrl T_s Período de muestreo 0,25/0,5 µs F_s Frecuencia de muestreo (Fs) 4/2 MHz t_sim Tiempo final de simulación 10 s plane_vel Velocidad de la aeronave 0,0694/0,25 km/s signal_vel Velocidad de propagación de la señal 2,9979e5 km/s Channel FSPLEnable Control atenuación manual/automático 0 (manual) ∗ AF1Attenuation factor del slave X 0 dB AF2Attenuation factor del slave Y 0 dB SNR Signal-to-noise ratio [8···40] dB Rx SIM_INI Inicio del funcionamiento del Rx 801 muestras AdapPwrEnable Control activación autoadaptación umbrales 1 (ON) ∗ PostProcEnable Control activación posprocesado 1 (ON) ∗ BW Ancho de banda 3,1831e3 Hz rho ρen (3.13) 1,4 ∗ eta ηen (3.14) 1,45 ∗ τ1Constante de tiempo lenta en DdE 50e-6 s τ2Constante de tiempo rápida en DdE 0,5e-6 s Path tray_type Control del tipo de trayectoria 1/2/3 ∗ Información de la Tabla 4.2 Loran_Rx W_pulse Ancho de envolvente de referencia 801 muestras W_1 Ventana de búsqueda de máximos 1000·Fsmuestras TH_MASTER_WIN_MAX Tamaño máx. de separación de último pulso master 2250·Fsmuestras TH_MASTER_WIN_MIN Tamaño mín. de separación de último pulso master 1750·Fsmuestras TH_SLAVE_WIN_MIN Tamaño mín. de separación de último pulso slave 750·Fsmuestras TH_SLAVE_WIN_MAX Tamaño máx. de separación de último pulso slave 1250·Fsmuestras RAW_DETECTION_WINDOW Ventana de actividad del MD 10·W_1 muestras INACTIVE_TH Ventana de medida del ruido durante ME b2000·Fscmuestras GRI_WIN Tamaño máximo de la diferencia TOAS−TOAM(GRI −900)·Fsmuestras DelayMaster Retraso aplicable a TOAM9000·Fsmuestras DelaySlave Retraso aplicable a TOAS7000·Fsmuestras CalibratedDelay Retraso aplicable por desfases introducidos por filtros 41 muestras Tabla 4.1: Parámetros de configuración establecidos para las simulaciones. Estación Coordenadas [km] Parámetros [µs] CD BTT ED GRI Master [0 0] ∗ ∗ ∗ 49000 Slave X [600 0] 13000 2,0014e3 1,5001e4 ∗ Y [0 1200] 29000 4,0028e3 3,3003e4 ∗ Tabla 4.2: Información del chain empleado. Los datos no corresponden a ningún chain que figure en las tablas del Anexo A. conceptuales entre exactitud y precisión. Las medidas se dice que son exactas si la agrupación se encuentra cercana al valor real, y se dice que son precisas si entre ellas están cerca. Por ejemplo, en la Fig. 4.2, el diagrama de arriba a la derecha muestra medidas inexactas pero precisas, mientras que en el de abajo a la izquierda muestra medidas exactas pero imprecisas. Estos conceptos es necesario tenerlos comprendidos cuando se trabaja con cualquier sistema de medida, pero en especial con los sistemas de posicionamiento. Matemática y estadísticamente, tanto la exactitud como la precisión en 2D se medirán aquí 51 con 4.11 y 4.13 con 4.14. Ya con la intención de disipar cualquier recelo sobre lo aquí expuesto, los valores numéricos de los ajustes normales se recopilan en las Tablas 4.3 y 4.4. Las filas comparables se han sombreado con la misma tonalidad de gris. Las conclusiones son las siguientes: ↓SNR =⇒↑ ˆσ⇐⇒ dilatación de campana ⇐⇒ pérdida de precisión ˆσoriginales >ˆσdiezmadas =⇒mejoría en la precisión con el diezmado (tanto para ATE como para XTE) ˆµoriginales >ˆµdiezmadas =⇒mejoría en la exactitud con el diezmado (para XTE) Al igual que pasaría para el caso del posicionamiento estático con Ts, en esta ocasión puede constatarse que en términos generales la dependencia de TD con el tiempo no es realmente significativa puesto que para 250 y 900 km/s se obtienen valores prácticamente iguales. Se justifica de esta manera que en la práctica puede ser asumido que P(tX)≃P(tY)≃P(tY+τ)(§3.4.1.1). Póngase en relieve que se ha obtenido una apreciable mejoría en la estimación mediante el uso de un simple filtro FIR de diezmado fácilmente implementable; sin necesidad de poner en juego ningún tipo de ecualizador. SNR [dB] 40 28 18 11 ˆµ0,03458 0,094707 0,31179 0,68955 XTE ˆσ0,037083 0,077913 0,23082 0,5203 ˆµ-0,047458 -0,044838 -0,031233 0,062144 Estimaciones originales ATE ˆσ0,057105 0,18315 0,64023 1,3218 ˆµ0,024269 0,035865 0,10332 0,14701 XTE ˆσ0,022353 0,045989 0,10834 0,12187 ˆµ-0,049165 -0,054729 -0,042263 0,044581 Estimaciones diezmadas ATE ˆσ0,029232 0,050175 0,14072 0,28542 Tabla 4.3: Valores estimados en km de la media y desviación estándar para el ajuste gaussiano de los errores en el caso de vuelo a 250 km/h. SNR [dB] 40 28 18 11 ˆµ0,034353 0,10237 0,30351 0,7321 XTE ˆσ0,032856 0,082835 0,22998 0,56553 ˆµ-0,053124 -0,031013 0,0086978 -0,052624 Estimaciones originales ATE ˆσ0,048033 0,18806 0,57983 1,4019 ˆµ0,016182 0,031258 0,061835 0,12987 XTE ˆσ0,011344 0,023755 0,042521 0,08263 ˆµ-0,050083 -0,027207 -0,05794 -0,12659 Estimaciones diezmadas ATE ˆσ0,016331 0,044342 0,16397 0,42969 Tabla 4.4: Valores estimados en km de la media y desviación estándar para el ajuste gaussiano de los errores en el caso de vuelo a 900 km/h. 58 295 296 297 298 299 300 301 302 303 304 305 Coordenada x [km] 598.5 599 599.5 600 600.5 601 601.5 602 602.5 603 Coordenada y [km] SNR = 11 dB SNR = 18 dB SNR = 28 dB SNR = 40 dB Tray. real Figura 4.7: Posiciones estimadas para diferentes valores de SNR volando a 250 km/h en trayectoria recta. 295 296 297 298 299 300 301 302 303 304 305 Coordenada x [km] 598 599 600 601 602 603 604 Coordenada y [km] SNR = 11 dB SNR = 18 dB SNR = 28 dB SNR = 40 dB Tray. real Figura 4.8: Posiciones estimadas para diferentes valores de SNR volando a 900 km/h en trayectoria recta. 59 300 300.5 301 599.9 600 600.1 600.2 600.3 600.4 600.5 Coordenada y [km] SNR = 40 dB Est. originales Est. diezmadas Tray. real 299.5 300 300.5 301 599.6 599.8 600 600.2 600.4 600.6 SNR = 28 dB Est. originales Est. diezmadas Tray. real 299 300 301 302 Coordenada x [km] 599 599.5 600 600.5 601 Coordenada y [km] SNR = 18 dB Est. originales Est. diezmadas Tray. real 295 300 305 Coordenada x [km] 598 599 600 601 602 SNR = 11 dB Est. originales Est. diezmadas Tray. real Figura 4.9: Posiciones estimadas diezmadas un factor 12 frente a las originales para diferentes valores de SNR volando a 250 km/h en trayectoria recta. 60 SNR = 40 dB 0 0.05 0.1 0.15 0 50 100 Frecuencia SNR = 40 dB -0.2 -0.1 0 0.1 0.2 0 20 SNR = 28 dB 0 0.1 0.2 0.3 0 50 Frecuencia SNR = 28 dB -0.5 0 0.5 0 20 40 SNR = 18 dB 0 0.5 1 0 20 40 Frecuencia SNR = 18 dB -1 0 1 0 20 SNR = 11 dB 0 1 2 XTE [km] 0 20 Frecuencia SNR = 11 dB -4 -2 0 2 4 ATE [km] 0 10 20 30 Figura 4.10: Histogramas y ajuste a distribución normal de errores ATE y XTE para diferentes valores de SNR volando a 250 km/h en trayectoria recta; estimaciones originales. 61 SNR = 40 dB 0 0.02 0.04 0.06 0.08 0 5 10 Frecuencia SNR = 40 dB -0.1 -0.05 0 0 5 10 SNR = 28 dB 0 0.05 0.1 0.15 0 5 10 Frecuencia SNR = 28 dB -0.2 -0.1 0 0 5 SNR = 18 dB 0 0.1 0.2 0.3 0.4 0 5 Frecuencia SNR = 18 dB -0.4 -0.2 0 0.2 0 5 SNR = 11 dB 0 0.2 0.4 XTE [km] 0 5 Frecuencia SNR = 11 dB -0.5 0 0.5 ATE [km] 0 5 Figura 4.11: Histogramas y ajuste a distribución normal de errores ATE y XTE para diferentes valores de SNR volando a 250 km/h en trayectoria recta; estimaciones diezmadas. 62 300 301 302 599.5 600 600.5 601 601.5 Coordenada y [km] SNR = 40 dB Est. originales Est. diezmadas Tray. real 300 301 302 303 599.5 600 600.5 601 601.5 602 SNR = 28 dB Est. originales Est. diezmadas Tray. real 298 300 302 304 Coordenada x [km] 599.5 600 600.5 601 601.5 602 Coordenada y [km] SNR = 18 dB Est. originales Est. diezmadas Tray. real 296 298 300 302 304 Coordenada x [km] 598 599 600 601 602 603 SNR = 11 dB Est. originales Est. diezmadas Tray. real Figura 4.12: Posiciones estimadas diezmadas un factor 12 frente a las originales para diferentes valores de SNR volando a 900 km/h en trayectoria recta. 63 SNR = 40 dB 0 0.05 0.1 0.15 0 50 100 Frecuencia SNR = 40 dB -0.2 -0.1 0 0.1 0 20 SNR = 28 dB 0 0.1 0.2 0.3 0 50 Frecuencia SNR = 28 dB -0.5 0 0.5 0 20 40 SNR = 18 dB 0 0.5 1 0 20 40 Frecuencia SNR = 18 dB -1 0 1 0 20 40 SNR = 11 dB 0 1 2 3 XTE [km] 0 20 40 Frecuencia SNR = 11 dB -4 -2 0 2 4 ATE [km] 0 20 40 Figura 4.13: Histogramas y ajuste a distribución normal de errores ATE y XTE para diferentes valores de SNR volando a 900 km/h en trayectoria recta; estimaciones originales. 64 SNR = 40 dB 0 0.02 0.04 0 5 Frecuencia SNR = 40 dB -0.08 -0.06 -0.04 -0.02 0 0 5 SNR = 28 dB 0 0.05 0.1 0 5 Frecuencia SNR = 28 dB -0.1 0 0.1 0 5 SNR = 18 dB 0 0.05 0.1 0.15 0 5 Frecuencia SNR = 18 dB -0.4 -0.2 0 0.2 0.4 0 5 SNR = 11 dB 0 0.1 0.2 0.3 XTE [km] 0 5 Frecuencia SNR = 11 dB -1 -0.5 0 0.5 1 ATE [km] 0 5 Figura 4.14: Histogramas y ajuste a distribución normal de errores ATE y XTE para diferentes valores de SNR volando a 900 km/h en trayectoria recta; estimaciones diezmadas. 65 4.2.2.2. Trayectoria circular Otro tipo de trayectorias efectuadas por las aeronaves con asiduidad son las trayectorias circulares; que en general suelen acotarse en arcos de circunferencia, llamados comúnmente virajes, y que sirven habitualmente para empatar un rumbo con otro. Se ha tomado un arco de circunferencia de radio R= 0,85 km, punto de origen (300,600) y con sentido horario. Obsérvese que en la Fig. 4.5 se han unido los puntos PkyPk+1 mediante una línea recta. Esto hace que para la trayectoria recta, los parámetros ATE y XTE funcionen adecuadamente a la hora de caracterizar el error. Sin embargo, adviértase que si la trayectoria es curva, los puntos PkyPk+1 no se encuentran unidos mediante una línea recta, sino por una curva, lo que induce una cierta inexactitud a la hora de caracterizar el error con estos parámetros. Sin embargo, la alta tasa de adquisición de estimaciones junto con el amplio radio con el que se suelen ejecutar los virajes hace que se puedan emplear ATE y XTE con sobrada tranquilidad, ya que la aproximación recta de esos puntos se aproxima mucho a la curva. Las estimaciones obtenidas para este escenario se hallan en la Fig. 4.15. Además se recopilan en la Tabla 4.5, igual que anteriormente, los valores numéricos de los ajustes tanto para las estimaciones originales como para las diezmadas. No es difícil comprobar que los valores numéricos de la tabla prácticamente coinciden con los mostrados para el caso de la trayectoria recta, lo que implica que la operatividad del sistema en ambos casos no difiere sustancialmente. SNR [dB] 40 28 18 11 ˆµ0,041561 0,091962 0,27979 0,69006 XTE ˆσ0,024933 0,069295 0,20067 0,53014 ˆµ-0,036742 -0,033479 -0,028708 0,13426 Estimaciones originales ATE ˆσ0,0598 0,18452 0,55166 1,3216 ˆµ0,035998 0,030174 0,063069 0,18573 XTE ˆσ0,017683 0,026084 0,057378 0,27322 ˆµ-0,035516 -0,025871 -0,036132 0,18017 Estimaciones diezmadas ATE ˆσ0,03447 0,066489 0,16645 0,34782 Tabla 4.5: Valores estimados en km de la media y desviación estándar para el ajuste gaussiano de los errores en el caso de vuelo a 900 km/h. 4.3. Efecto del canal 4.3.1. Ruido El canal AWGN se comporta de manera tal que su salida yen un tiempo discreto kes la suma de la entrada xen ese mismo instante más el ruido z: y[k] = x[k] + z[k](4.5) zes una variable aleatoria distribuida acorde a una normal de media cero y varianza σ2 n{z[k]∼ N(0,σ2 n)}. Además, z[k]tiene la propiedad de ser independiente (autocorrelación nula) y estar idénticamente distribuido respecto al ruido de instantes 66 296 297 298 299 300 301 302 303 304 305 Coordenada x [km] 599 599.5 600 600.5 601 601.5 602 602.5 603 Coordenada y [km] SNR = 11 dB SNR = 18 dB SNR = 28 dB SNR = 40 dB Tray. real Figura 4.15: Posiciones estimadas para diferentes valores de SNR volando a 900 km/h en viraje. diferentes. El comando randn de matlab devuelve valores escalares aleatorios ajustados a la distribución normal estándar [20]. Por lo tanto: u←randn |u∼ N(0,1). La propiedad de simetría respecto a su media de la distribución normal permite relacionar cualquier variable aleatoria normal v∼ N(µ,σ2)con la distribución normal estándar ude la siguiente manera: v=u·σ+µ(4.6) Si particularizamos para nuestra variable aleatoria z=v, se obtiene fácilmente que z=u·σn(4.7) La cuestión ahora es qué valor debe de adoptar σnpara que la señal recibida a la entrada del transmisor presente un valor de SNR fijado. Para ello se comienza con la definición del SNR: snr =Px Pn,Rx (4.8) donde Pxes la potencia de la señal transmitida y Pn,Rx es la potencia del ruido a la entrada del receptor. SNR es simplemente snr en dB. La potencia de una señal en un determinado intervalo de frecuencias [f0,f1]puede ser calculado integrando la densidad espectral de potencia Sxx(f)con respecto a la frecuencia. Como Sxx(f)es simétrica respecto a f= 0 la misma cantidad de potencia puede ser atribuida a frecuencias tanto positivas como negativas, lo que provoca la aparición de un factor 2 [18]. P= 2 Zf1 f0 Sxx(f)df (4.9) 67 centésimas de kilómetro, o lo que es lo mismo, decenas de metros. Cuando se aumenta progresivamente la hostilidad del canal empeoran los errores cometidos, llegando incluso a varios kilómetros. Con la finalidad de mejorar los datos obtenidos se ha empleado un proceso de diezmado mediante un filtro FIR para el caso de posicionamiento móvil, con el que se ha logrado reducir sustancialmente los errores. En concreto, se consigue disminuir la dispersión de las estimaciones en el entorno del 70 %. 5.2. Propuestas de ampliación A pesar de haberse cumplido el objetivo principal propuesto para este trabajo, siempre existen aspectos que son susceptibles de ampliación y mejora. A continuación se reúnen aquellos que han surgido a lo largo del desarrollo del mismo y que se creen interesantes de abordar en un futuro con el objetivo de enriquecer el estudio aquí comenzado. El que sin duda sería el aspecto más interesante a modelar sería la inclusión de chains que dispusieran de más de dos estaciones slaves. Esto haría posible que, en principio y a expensas del estudio aparejado, se obtuviera una mejora en la estimación de la posición ya que esta se pasaría a obtener a partir de más de una única pareja de TD como se ha hecho hasta ahora. La redundancia de información a la hora de ofrecer una estimación favorece por lo general la fiabilidad del resultado alcanzado. Para su resolución sería necesario la generalización del método de Chan para el caso de M≥4, extensión que se encuentra incluida en [21]. A pesar de que esto infiera un aumento en la complejidad matemática involucrada, la programación los códigos de matlab se ha orientado de manera que no resulte un arduo trabajo la extensión a un mayor números de estaciones slaves, incitando a abordar esta sugerencia. Incluso la integración de este sistema con su candidato idóneo, GPS, por medio de un ecualizador del estilo del filtro de Kalman satisfaría la tendencia habitual de integrar sistemas de navegación con el fin de mejorar la precisión de las estimaciones. Nótese que se ha modelado un único chain y por lo tanto, una cuestión importante que ha quedado sin resolver: ¿cómo afecta el CRI que surge cuando las áreas de cobertura de dos o más chains se solapan? Este tipo de interferencia derivado del uso de la misma frecuencia para todas las estaciones motivó el surgimiento de la diferencia en valor del GRI para cada uno de los chains con propósito de poder identificarlos. Poner en juego más de un chain y ser capaz de identificarlos y distinguirlos conlleva perspicaces actuaciones en el bloque del receptor. Adicionalmente con esta disposición de más de un chain, podría también modificar el receptor propuesto para que la sintonización del chain que proporcione la mejor estima de la posición se realizara automáticamente, evitando así solicitar la actuación por parte del navegante. A pesar de que en este trabajo se haya optado por utilizar el método de Chan como estimador hiperbólico, nada impide que pueda emplearse otro que incluso pueda ser mejor. La coexistencia de varios estimadores hiperbólicos en la literatura no suele ser caprichosa, sino que atiende a que cada uno de ellos se comportará de una manera específica en función del escenario concreto en 74 el que operen. Se propone el uso de cualesquiera otros y su comparativa con el aquí propuesto para ilustrar sus diferencias. Los efectos en los cambios en los medios de propagación, así como la presencia de las sky waves y las correspondientes acciones correctivas, podrían incluirse con la intención de reproducir el comportamiento real del sistema con mayor rigor. Igualmente, el tratamiento del ruido introducido por el canal como un AWGN podría extenderse para introducir otros tipos de ruido (CRI y/o CRW por ejemplo). 75 Nomenclatura HMSemiplano master HSSemiplano slave AEE Average Euclidean Error AF Attenuation Factor AFE Analogue Front End ASF Additional Secondary phase Factor ATE Along-Track Error AWGN Additive White Gaussian Noise BPF Band-Pass Filter BTT Baseline Travel Time BW Bandwidth CC Cross-Correlation CD Coding Delay CEP Circular Error Probability CRI Cross-Rate Interference CWI Carrier-Wave Interference DAC Divide And Conquer DdE Detector de Envolvente DRMS Distance Root Mean Square ED Emission Delay eLORAN Enhanced LORAN FAA Federal Aviation Administration FGDC Federal Geographic Data Committee FIR Finite Impulse Response 77 FSPL Free-Space Path Loss GNSS Global Navigation Satellite System GPS Global Positioning System GRI Group Repetition Interval HEA Harbor Entrance and Approach LDC Loran Data Channel LF Low Frequency LOP Line Of Position LORAN Long Range Navigation LPF Low-Pass Filter MD Modo Detección ME Modo Escucha NPA Non-Precision Approach PC Phase Coding PCI Phase-Code Interval PE Position Estimator PF Primary phase Factor PGEN Pulse Generator PRF Pulse Repetition Frequency PRR Pulse Repetition Rate PRT Position Refresh Time PSD Power Spectral Density RC Receiver Control Rx Receptor SD Signal Detector SF Secondary phase Factor SI Spherical-interpolation SNR Signal-to-Noise Ratio SX Spherical-intersection 78 TAC Transmitter Automatic Controller TD Time Difference TOA Time Of Arrival TTG Time To Go Tx Transmisor USCG United States Coast Guard XTE Cross-Track Error 79 Bibliografía [1] L.A. Turner, J.A. Pierce, A.A. McKenzie & R.H. Woodward. LORAN. 1.aed. McGraw Hill, 1948. [2] G.W. Johnson, P.F. Swaszek, R.J. Hartnett, R. Shalaev & M. Wiggins. «An Evaluation of eLoran as a Backup to GPS». En: 2007 IEEE Conference on Technologies for Homeland Security. Mayo de 2007, págs. 95-100. doi:10. 1109/THS.2007.370027. [3] John A. Volpe National Transportation Systems Center. Vulnerability Assessment of the Transportation Infrastructure relying on the Global Positioning System. Inf. téc. Autor, agosto de 2001. [4] International Loran Association. Enhanced Loran - Definition Document. Octubre de 2007. [5] M. Narins. «The Global Loran / eLoran Infrastructure Evolution». En: 2014 National Space-Based Positioning, Navigation, and Timing Meeting. Junio de 2014. [6] U.S. Coast Guard. «Terminate Long Range Aids to Navigation (Loran-C) Signal». En: Federal Register 75.4 (Enero de 2010). Docket NoUSCG-2009–0299, pág. 998. [7] Durk Van Willigen. «Eurofix». En: Journal of Navigation 42.3 (1989), págs. 375-381. doi:10.1017/S0373463300014661. [8] L. Tetley & D. Calcutt. «Electronic Navigation Systems». En: 3.aed. Butterworth- Heinemann, 2001. Cap. 4. [9] J. de Burgos. Álgebra lineal y geometría cartesiana. 3.aed. Madrid: McGraw- Hill España, 2013. [10] National Geospatial-Intelligence Agency. «Loran Navigation». En: American Practical Navigator: An epitome of navigation. Vol. 1. Springfield, Virginia, EE.UU.: Autor, 2017. Cap. 24. [11] Navy Department. The Tactical Use of Radar in Aircraft. Radar Bulletin 2A. Washington: Office of the Chief of Naval Operations, 1946. [12] U.S. Coast Guard. «Loran System». En: Aids to Navigation manual. Washington, EE.UU.: Autor, enero de 1953. Cap. 18. [13] U.S. Coast Guard. Specification of the transmitted Loran-C signal. Commandant Instruction M16562.4A. Washington D.C., U.S.A.: U.S. Dept. of Transportation, mayo de 1994. 81 [14] R.L. Frank. «Multiple Pulse and Phase Code Modulation in the Loran-C System». En: IRE Transactions on Aeronautical and Navigational Electronics 7.2 (Junio de 1960), págs. 55-61. issn: 0096-1639. doi:10.1109/TANE3.1960. 4201731. [15] A.J. Fisher. The Loran-C Cycle Identification Problem. Department of Computer Science, University of York, U.K., 1999. [16] J.L. Besada & M. Sierra. Tema 3: Propagación de ondas en medio natural. OpenCourseWare de la Universidad Politécnica de Madrid: Radiación y Propagación. Madrid, 2010. [17] H.T. Sherman & V.L. Johnson. The Loran-C Navigation Station. Inf. téc. Institute of Navigation, National Marine Meeting, octubre de 1975. [18] Wikipedia, la enciclopedia libre.url:https://es.wikipedia.org/wiki/. [19] C.W. Dull & W.H. Clerke. «Loran-C Navigation as an Aid to Aerial Photogaphic Operations». En: Photogrammetric Engineering and Remote Sensing 49.11 (Noviembre de 1983), págs. 1577-1584. [20] MathWorks c .Matlab Documentation.url:https://es.mathworks.com/ help/matlab/. [21] Y.T. Chan & K.C. Ho. «A simple and efficient estimator for hyperbolic location». En: IEEE Transactions on Signal Processing 42.8 (Agosto de 1994), págs. 1905-1915. issn: 1053-587X. doi:10.1109/78.301830. [22] G.A. Mizusawa. «Performance of hyperbolic position location techniques for code division multiple access». Tesis de mtría. Blacksburg, Virgina, EE.UU.: Virginia Polytechnic Institute and State University, agosto de 1996. [23] M. Beckmann. «Carrier Wave Signals Interfering with Loran-C». Tesis doct. Technische Universiteit Delft, enero de 1992. [24] L. Lee, M. Jones, G.S. Ridenour, M.P. Testa & M.J. Wilson. «Investigating and Comparing Spatial Accuracy and Precision of GPS-Enabled Devices in Middle Tennessee». En: Geo-Informatics in Resource Management and Sustainable Ecosystem: Second International Conference, GRMSE 2014. Ed. por F. Bian & Y. Xie. Springer-Verlag Berlin Heidelberg, 2015, págs. 215-224. doi:10. 1007/978-3-662-45737-5_22. [25] GPS World staff. GPS Accuracy: Lies, Damn Lies, and Statistics. Enero de 1998. url:http://www.gpsworld.com/gps-accuracy-lies-damn-lies- and-statistics/. [26] Federal Aviation Administration. RNP 10 - Qualification through Data Collection. Autor. Abril de 2019. [27] D.A. Panshin. «What You Should Know About Loran-C Receivers». En: Marine Electronics. Washington D.C., U.S.A.: ERIC Clearinghouse, Mayo de 1978. 82 Anexos 83 Anexo B Método de Chan B.1. Modelo matemático general para la estimación hiperbólica de la posición Se considera el plano bidimensional cartesiano en el que se encuentran situados tanto el vehículo del que se desea estimar su posición, con coordenadas (x,y), como las diferentes estaciones emisoras que conforman un chain del sistema Loran-C (tanto slaves como master), cuyas coordenadas son (xj,yj). Por tanto, la distancia entre el vehículo y cualquier estación del chain en un instante determinado será rj=q(xj−x)2+ (yj−y)2para j=M,X,Y,Z... (B.1) Ya que todas las estimaciones del TDi(i=X,Y,Z...)ofrecidas por el receptor se encuentran referidas a la estación master, la diferencia entre la distancia del vehículo a la estación i-ésima y a la estación master es ri,M := ri−rM=     −2asi ri< rM(semiplano slave) 0si ri=rM(baseline bisector) 2asi ri> rM(semiplano master) (B.2) donde hay que recordar que a=a(TD), por lo que ri,M =ri,M (TD), siendo por tanto esta una medida (dato) proporcionada por el receptor. Las ecs. (B.2) definen un conjunto de ecuaciones hiperbólicas no lineales cuya solución (intersección entre ellas) proporciona las coordenadas del vehículo. Sin embargo, la resolución de estas ecuaciones no lineales es compleja. Por ello, una posible (y usual) opción es linealizar (B.2) mediante expansiones por series de Taylor reteniendo los dos primeros términos y posteriormente resolver iterativamente hasta converger a la solución (Método de Foy, [22]). Un método alternativo frecuentemente utilizado consiste en transformar este mismo conjunto de ecuaciones no lineales en otro conjunto de ecuaciones. Este segundo enfoque será el llevado a cabo. Por tanto, despejando ride (B.2) y elevando al cuadrado se tiene r2 i= (ri,M +rM)2(B.3) en la que se introduce (B.1) y se desarrollan términos, obteniendo r2 i= (xi−x)2+ (yi−y)2= (ri,M +rM)2(B.4) Ki+x2−2xix+y2−2yiy=r2 i,M +r2 M+ 2ri,M rM(B.5) 91 donde Ki:= x2 i+y2 i. Si ahora le restamos r2 Ma (B.5) se tiene Ki+x2−2xix+y2−2yiy−r2 M=r2 i,M + 2ri,M rM(B.6) e introducimos la ec. (B.1) particularizada para j=M, la ecuación resultante sería r2 i,M + 2ri,M rM=Ki−KM−2xxi,M −2yyi,M con i=X,Y,Z... (B.7) en la que xi,M := xi−xMyyi,M := yi−yM. Nótese que las únicas incógnitas que aparecen son x,yyrM, ya que: Sabiendo las coordenadas (xi,yi)de todas las estaciones, se pueden conocer las cantidades Kj,xi,M yyi,M . El receptor, a través de la estimación del TD, proporcionará ri,M . El conjunto de ecs. (B.7) son ahora lineales, lo que facilita en gran medida su resolución. B.2. Propuesta de resolución del Método de Chan Para el caso concreto de tener una estación master (M) y dos estaciones slaves (SXySY), el problema tiene solución analítica, ya que tenemos el mismo número de incógnitas que de ecuaciones ((B.7) junto con (B.1)). El Método de Chan parte de (B.7), que particularizada para esta composición se tiene i=X=⇒r2 X,M + 2rX,M rM=KX−KM−2xxX,M −2yyX,M (B.8) i=Y=⇒r2 Y,M + 2rY,M rM=KY−KM−2xxY,M −2yyY,M (B.9) Si se despeja xey(vector de posición) en función de rMy se expresa en forma matricial, el resultado quedaría ~ x(rM) = x(rM) y(rM)=ArM+B(B.10) donde M=−xX,M yX,M xY,M yY,M −1 (B.11) A=MrX,M rY,M (B.12) B=1 2Mr2 X,M −KX+KM r2 Y,M −KY+KM(B.13) Para garantizar que el sistema tiene solución única (sistema compatible determinado) se asume la existencia de la inversa de la matriz M, es decir, tiene rango máximo. Esto se cumple siempre que las estaciones no se encuentren emplazadas linealmente [21]. Si se despeja cada componente de (B.10) x y=a1 a2rM+b1 b2=⇒(x=a1rM+b1 y=a2rM+b2 (B.14) 92 y se introduce en (B.1) con j=M rM=p(xM−(a1rM+b1))2+ (yM−(a2rM+b2))2(B.15) r2 M=KM+(a1rM+b1)2−2xM(a1rM+b1)+ (a2rM+b2)2−2yM(a2rM+b2)(B.16) Si se reúnen términos en torno a la variable independiente rMen la ec. (B.16) se puede reconocer fácilmente que se trata de una simple ecuación cuadrática αr2 M+ βrM+γ= 0, donde los coeficientes valen α=a2 1+a2 2−1(B.17) β= 2a1b1−2xMa1+ 2a2b2−2yMa2(B.18) γ=KM−2xMb1+b2 2−2yMb2+b2 1(B.19) De las dos posibles raíces que puede tener a lo sumo la ecuación cuadrática anterior, sólo tendrán sentido físico aquellas que cumplan que rM={rM∈R:rM≥0}, que introducidas de vuelta en (B.10) facilitan finalmente la estimación de la posición ~x. En caso de que ambas raíces cumplieran dicha condición, se encontraría ante una ambigüedad que habría que resolver estableciendo algún criterio que nos indicara cuál de las dos soluciones es la más probable a ser la solución real. En esta ocasión, este criterio de selección se ha establecido para que se tome como estimación acertada aquella que diste menos de la estimación anteriormente obtenida. Si ocurriese que la primera de las estimaciones ya es ambigua, al no disponer de una ‘estimación anterior’, se considerará que será conocida una estimación preliminar aproximada (obtenida por algún medio externo al sistema Loran-C) al inicio del funcionamiento del receptor con la que poder cotejarla. Sin embargo, en caso de que el chain disponga de 3 o más estaciones slaves, el sistema (B.7) se encontrará sobredeterminado, ya que se dispondrán de un mayor número de ecuaciones que de incógnitas. En la práctica, y de manera acentuada con la inclusión de ruido, es improbable que el sistema sea compatible. Por ello, la respuesta apropiada a la búsqueda de la estimación de la posición será aquella pareja (x,y)que mejor se ajuste, por ejemplo en el sentido de los mínimos cuadrados, a las ecuaciones. 93