1
METODOLOGÍA PARA LA REALIZACIÓN EN PARALELO DE FILTROS DE ORDEN SUPERIOR SOBRE PROCESADORES DIGITALES DE SEÑALES
ANDRÉS FELIPE JARAMILLO A.
UNIVERSIDAD TECNOLÓGICA DE PEREIRA FACULTAD DE INGENIERÍAS
PROGRAMA DE INGENIERÍA FÍSICA PEREIRA, COLOMBIA
2015
2
METODOLOGÍA PARA LA REALIZACIÓN EN PARALELO DE FILTROS DE ORDEN SUPERIOR SOBRE PROCESADORES DIGITALES DE SEÑALES
ANDRÉS FELIPE JARAMILLO A.
Proyecto de grado para optar por el título de Ingeniero Físico
DIRECTOR: M.Sc EDWIN ANDRES QUINTERO SALAZAR
DOCENTE DEL PROGRAMA DE INGENIERÍA FISICA
UNIVERSIDAD TECNOLÓGICA DE PEREIRA FACULTAD DE INGENIERÍAS
PROGRAMA DE INGENIERÍA FÍSICA PEREIRA
2015
3
Nota de aceptación: ___________________________________ ___________________________________ ___________________________________ ___________________________________ ___________________________________ Firma del presidente del jurado ___________________________________ Firma del jurado Pereira, 2015
4
AGRADECIMIENTOS A mi padre Leonardo Jaramillo Z. por inculcarme el amor a la ciencia, a mi madre Miriam Eugenia Alvarado R. por enseñarme que para salir adelante se necesita firmeza y convicción, a mi hermana Adriana María Santa A. por regalarme el valor de la responsabilidad y a mi hermano Ricardo Andrés Jaramillo A. por no dejarme olvidar el niño que llevo dentro. A mis amigos Jorge Alejandro Bueno V., Miguel Ángel Jaramillo S. por estar conmigo en todo momento, A mi novia Caterine Tatiana Cifuentes R. por acompañarme en mi camino y ser una fuente de apoyo incondicional para mi progreso. Un agradecimiento especial para el Ingeniero Edwin Andrés Quintero quien ha formado parte fundamental de mi formación académica y personal. A mi familia, a mis amigos, a mis profesores, a mis compañeros de carrera y a la vida misma por hacerme quien soy.
5
ÍNDICE GENERAL
Contenido ÍNDICE DE FIGURAS ...................................................................................................................... 7
ÍNDICE DE TABLAS ........................................................................................................................ 8
RESUMEN ......................................................................................................................................... 9
INTRODUCCIÓN ............................................................................................................................ 10
1. DESCRIPCIÓN DE LA PROPUESTA ................................................................................ 11
1.1. Definición del Problema ............................................................................................. 11
1.2. Justificación .................................................................................................................. 12
1.3. Objetivos ......................................................................................................................... 13
1.3.1. Objetivo General ........................................................................................................... 13
1.3.2. Objetivos Específicos ................................................................................................. 13
1.4. Estado del Arte.............................................................................................................. 14
2. SEÑALES Y SISTEMAS ....................................................................................................... 18
2.1 Señales Analógicas ......................................................................................................... 18
2.2 Señales Discretas ............................................................................................................. 18
2.3 Sistemas Físicos .............................................................................................................. 19
2.4 Sistemas en tiempo continuo ....................................................................................... 19
2.5 Sistemas en tiempo discreto ......................................................................................... 20
2.6 Relación entre sistemas en tiempo continuo y discreto ....................................... 21
2.7 Bancos de filtros [15] ...................................................................................................... 22
3. PROCESADORES DIGITALES DE SEÑALES Y HERRAMIENTAS
ARQUITECTÓNICAS .................................................................................................................... 23
3.1 Procesadores digitales de señales (DSPs) ............................................................... 23
3.2 ADSP Blackfin-533 [16] ................................................................................................... 23
3.3 Arquitectura de BUS y Memoria de la ADSP BF-533 .............................................. 25
3.4 Archivos de registros y unidades de procesamiento de Hardware .................... 27
3.5 Periféricos Básicos. ......................................................................................................... 28
3.6 Entorno de Desarrollo Integrado. ................................................................................ 29
3.7 Visual DSP ++ 5.0 [18] .................................................................................................... 29
3.8 Manejo de proyectos sobre Visual DSP ++ 5.0 ......................................................... 30
6
3.9 Depuración sobre Visual DSP ++ 5.0 .......................................................................... 32
3.10 Herramientas de medición: Tarjetas de Adquisición de datos. ........................... 34
3.11 NI MyDAQ [19]. .................................................................................................................. 35
3.12 Estructura física de la NI MyDAQ. ................................................................................ 35
3.13 Componentes del Hardware de la NI MyDAQ. .......................................................... 36
3.14 Funcionalidades Pre-programadas para la NI MyDAQ. ......................................... 38
4. METODOLOGÍA ..................................................................................................................... 42
4.1 Perspectiva General y Consideraciones Iniciales ................................................... 42
4.2 Base Teórica ...................................................................................................................... 42
4.3 Definiciones y Convenciones ....................................................................................... 43
4.4 Algoritmos y Flujos.......................................................................................................... 44
4.5 Implementación sobre Visual DSP++ .......................................................................... 48
4.6 Consecuencias Prácticas ............................................................................................... 56
5. RESULTADOS ....................................................................................................................... 58
6. CONCLUSIONES ................................................................................................................... 64
7. BIBLIOGRAFÍA ...................................................................................................................... 72
8. ANEXOS .................................................................................................................................. 75
7
ÍNDICE DE FIGURAS Figura 1. Hardware de la ADSP-BF533. ............................................................................................. 24
Figura 2. Modelo de memoria Jerárquico de la ADSP BF-533........................................................... 25
Figura 3. Arquitectura de memoria de la ADSP-BF533. .................................................................... 25
Figura 4. Mapa de la memoria de la ADSP-BF5xx. ............................................................................ 26
Figura 5. Núcleo de la ADSP-BF533 [17]. ........................................................................................... 27
Figura 6. Visual DSP++ 5.0. ................................................................................................................ 30
Figura 7. Ventana para la creación/edición de Sesiones en Visual DSP++. ....................................... 31
Figura 8. Cuadro para la configuración de nuevo proyecto en Visual DSP++. .................................. 32
Figura 9. Configuración para la depuración de un proyecto. ............................................................ 33
Figura 10. Procedimiento para la utilización de la herramienta de monitoreo gráfico de variables.
........................................................................................................................................................... 34
Figura 11. Tarjeta NI MyDAQ. ........................................................................................................... 35
Figura 12. Lateral de la NI MyDAQ. ................................................................................................... 35
Figura 13. Conectores de la NI MyDAQ para el uso del multímetro digital. (1) Medición de:
Resistencia, voltaje, diodo y continuidad. (2) Medición de corriente. ............................................. 36
Figura 14. Diagrama de bloques de los Componentes de la NI MyDAQ. .......................................... 36
Figura 15. Interfaz para el uso del multímetro digital de la NI MyDAQ. ........................................... 38
Figura 16. Interfaz para el uso de la funcionalidad: Osciloscopio digital de la NI MyDAQ. .............. 39
Figura 17. Interfaz para el uso de la funcionalidad: Generador de funciones de la NI MyDAQ. ...... 40
Figura 18. Interfaz para el uso de la funcionalidad: Analizador de Bode de la NI MyDAQ. .............. 41
Figura 19. Interfaz para el uso de la funcionalidad: Analizador dinámico de señales de la NI
MyDAQ. ............................................................................................................................................. 41
Figura 20. Algoritmo general de la implementación sobre la DSP. ................................................... 44
Figura 21. Estructura genérica del banco de filtros para la reconstrucción del sistema original. .... 57
8
ÍNDICE DE TABLAS Tabla 1. Transformaciones para llevar sistemas en tiempo continúo a tiempo discreto. ................ 21
Tabla 2. Herramientas para la depuración de código en un proyecto sobre Visual DSP++. ............. 33
Tabla 3. Cálculo del uso de memoria de la implementación por medio de la metodología
tradicional y la planteada. ................................................................................................................. 57
9
RESUMEN
Los filtros hacen parte integral de la gran mayoría de soluciones tecnológicas que actualmente se implementan, por ello, se hace importante realizar investigaciones en esta área respecto a la eficiencia de ellos. En este trabajo se presenta una metodología cuya implementación optimiza el uso del hardware disponible para el procesamiento. Esta consta de reconstruir el sistema de filtrado original mediante subsistemas los cuales son hallados mediante el algoritmo de fracciones parciales y, donde la superposición de las respuestas en frecuencia de estos, reconstruyen efectivamente el sistema original, además, todas las operaciones necesarias para el filtrado de la información se llevan a cabo mediante la optimización mencionada. La metodología fue implementada sobre la ADSP BF-533 donde se encontraron problemas para la realización de sistemas con polos anidados y donde los coeficientes que el numerador de la función de
transferencia fuesen del orden de ≈ 10−8. Se obtuvieron resultados satisfactorios para la reconstrucción de diferentes tipos de sistemas de filtrado de órdenes entre cuatro y diez con errores porcentuales en las frecuencias de corte medidas entre 0.23% y 11.9% con lo cual la metodología propuesta demostró su funcionalidad.
ABSTRACT
Generaly the filters are an integral part of technological solutions currently
being implemented, therefore it is important to do research in this area regarding the efficiency of them. This paper presents a methodology whose implementation optimizes the use of available hardware for processing is presented. This consists of reconstructing the original system filtrate by subsystems which are found by the algorithm of partial fractions and, where overlapping frequency responses of these, effectively reconstruct the original system, also all the operations necessary for filtering information are performed by said optimization. The impementacion of methodology was held on the ADSP BF-533 where found that implementation had problems with systems where has nested poles and the numerator of the transfer function have very too small coeficients. The results shown satisfactory for the reconstruction of several types of filter systems’s with orders from four to ten with percentage errors in the cutoff frequencies measured between 0.23% and 11.9% which the proposed methodology proved its functionality.
10
INTRODUCCIÓN
En la actualidad el ser humano se encuentra rodeado de innumerables implementaciones tecnológicas que realizan lo más eficientemente posible su trabajo. Las aplicaciones que llevan consigo sistemas de filtrado en su gran mayoría necesitan una alta selectividad frecuencial, esto conlleva ciertos problemas de implementación como lo son tiempos de procesamiento altos, uso de memoria elevado, restricciones de hardware, etc. Por esta razón día a día la investigación sobre la eficiencia de las implementaciones de sistemas de filtrado se hace más importante puesto que su influencia sobre la cotidianidad de las personas es cada vez más relevante por el gran flujo de aplicaciones que a diario se usan.
De allí que el desarrollo de una metodología para la implementación de
filtros que optimice el desempeño con respecto a las implementaciones tradicionales es importante para el futuro desarrollo de aplicaciones cada vez más robustas. En este trabajo se presenta una metodología de procesamiento con el propósito de brindar la posibilidad de optimizar el desempeño en la gran mayoría de aplicaciones de software o hardware que lleven consigo sistemas de filtrado LTI.
La metodología aquí presentada se implementó sobre la ADSP BF-533
donde mediante el uso de sus periféricos y el códec AD1836, se filtran las señales de entrada que llegan a la DSP por los puertos RCA.
11
1. DESCRIPCIÓN DE LA PROPUESTA
1.1. Definición Del Problema
En la instrumentación física los sistemas de filtrado juegan un rol fundamental
ya que estos cumplen con el objetivo de remover el ruido que posean las señales
provenientes de los sensores. En algunas aplicaciones se requiere un orden
elevado en el sistema de filtrado, ya que estas pueden requerir una alta
selectividad frecuencial y la forma comúnmente empleada para aumentar esta
característica es elevar el orden del sistema. Generalmente las implementaciones
de los sistemas de filtrado se hacen en el mundo digital puesto que las
realizaciones analógicas son tediosas y presentan variaciones por muchos
factores. Sin embargo, las realizaciones digitales presentan desventajas frente al
procesamiento analógico, de ellas, la latencia es la que presenta mayores
consecuencias. Ella consiste en un retardo notorio en el procesamiento de los
datos y a raíz de este fenómeno, se presenta pérdida de datos y por lo tanto,
errores al momento del filtrado de la señal. Esta problemática se presenta
principalmente cuando los sistemas de filtrado son robustos y de orden elevado.
Por otra parte en la actualidad el procesamiento concurrente [10] es una solución
viable cuando se requiere procesar una gran cantidad de datos de forma rápida y
eficiente. De allí y considerando las ventajas prácticas que ofrece el
procesamiento digital de señales y el aumento de eficiencia que proveen los
sistemas concurrentes, se nota la necesidad de responder la siguiente pregunta
de investigación:
- ¿Cómo construir e implementar una metodología que permita realizar filtros
de orden superior y optimizar su implementación sobre procesadores
digitales de señales?
12
1.2. Justificación
La necesidad actual de implementar sistemas de filtrado robustos y
complejos es producto del flujo en aumento de las aplicaciones sobre sistemas
embebidos. Por tal razón es necesario buscar la forma de disminuir la carga de
procesamiento sobre estos sistemas, lo cual se lograría si los filtros allí
implementados disminuyen su costo de procesamiento. En la actualidad el
procesamiento concurrente es una de las herramientas más utilizadas cuando se
requiere procesar una gran cantidad de datos de forma rápida y eficiente. Este tipo
de procesamiento se basa en una filosofía de paralelización de tareas. De allí se
notó la necesidad de aplicar esta filosofía sobre los sistemas digitales de filtrado
con la cual se obtendría una mejora en el rendimiento y se solventaría el problema
de la latencia para sistemas de filtrado robusto. Además, se podría brindar la
posibilidad al sistema que desarrolle su actividad en tiempo real. De esta forma se
espera implementar sistemas de filtrado muy robustos y complejos de forma
eficiente, sencilla y en tiempo real para el filtrado de señales sobre procesadores
digitales de señales.
13
1.3. Objetivos
1.3.1. Objetivo General
Diseñar e implementar una metodología para la realización en paralelo de
filtros lineales e invariantes en el tiempo que permita su ejecución en tiempo
real sobre procesadores digitales de señales.
1.3.2. Objetivos Específicos
Recopilar documentación e información sobre la arquitectura,
funcionamiento y tópicos prácticos sobre los procesadores digitales de
señales.
Desarrollar desde la parte teórica una metodología que permita
implementar filtros de orden superior en paralelo mediante ecuaciones en
diferencias.
Implementar la metodología desarrollada sobre el procesador digital de
señales Analog Device ADSP-BF533.
Validar los resultados esperados con mediciones adquiridas a partir de la
implementación de la metodología en el procesador digital de señales
Analog Device ADSP-BF533.
14
1.4. Estado del Arte
En la instrumentación física cualquiera que sea su aplicación, en general, se
tiene embebido una sección dedicada al procesamiento de las señales adquiridas
la cual permita el análisis posterior de dichas señales. Por esta razón se hace
fundamental que cada uno de los componentes del procesamiento sean
algoritmos, funciones, hardware etc. sean implementados de tal forma que su
eficiencia sea maximizada para así brindarle al sistema el mejor desempeño
posible. Para cumplir con este propósito a lo largo de la historia se han
desarrollado múltiples herramientas derivadas de diferentes investigaciones, las
cuales brindan al investigador la posibilidad de dotar a los sistemas de
procesamiento las características mencionadas. La invención de los Procesadores
Digitales de Señales (en adelante DSPs) revolucionó el mundo ya que permitió
tener en un sistema compacto todas las herramientas para realizar el
procesamiento de señales de una forma eficiente y robusta. Seguido se mostrarán
trabajos de otros investigadores relacionados a la base teórico/practica y el
enfoque que este proyecto presenta.
A mediados de la década de los 80’s las DSPs ya ocupaban lugares en
importantes aplicaciones por lo compacto de su implementación y las
características de desempeño que ofrecían. Una de las aplicaciones con mayor
interés para los investigadores de aquel entonces, era el reconocimiento de
palabras habladas por un usuario. En 1986 David O. Morgan y Kevin W. Leary, de
Brown University y de One Technology Way respectivamente, diseñaron una
metodología basada en el algoritmo DUR (Discrete Utterance Recognition) en [1],
con el fin de reconocer palabras habladas y la cual consistía de dos grandes
etapas. La primera constaba de diezmar los posibles candidatos mediante una
búsqueda rápida de las mejores coincidencias estadísticas. La segunda etapa,
mediante programación dinámica realizaba una búsqueda completa sobre los
candidatos restantes y elegía el candidato más acertado. Dicha metodología se
implementó sobre una ADSP-2100 de la Analog Device® la cual ofreció un
rendimiento adecuado para la aplicación. La metodología basada en dividir la
selección del candidato en dos etapas fue acertada, reduciendo el tiempo de
procesamiento y teniendo una alta tasa de certeza en el resultado. Los autores
manifestaron la necesidad de investigar otros enfoques basados en metodologías
de dos etapas para el reconocimiento de palabras puesto que aquí solo se exploró
una de ellas.
Desde la creación de las DSPs y hasta la actualidad uno de los problemas
centrales en su utilización ha sido el aumento del rendimiento de cada uno de los
componentes de la realización. Una de las posibles soluciones para aumentar el
15
rendimiento es la agrupación de varias DSPs trabajando en conjunto. En 1986
Katsuhiko Hayashi, Kaushal K. Dhar, Kazunori Sugahara y Kotaro Hirano
miembros de la IEEE propusieron varios esquemas para la realización de filtros
complejos mediante un grupo de DSPs interconectadas [2]. Plantearon la
descomposición de la función de transferencia H(z) en subsistemas e
implementarlos en paralelo sobre cada una las DSPs. Realizaron la comparativa
de rendimiento de la implementación de forma secuencial y en paralelo de H(z)
concluyendo que una implementación mediante N sistemas de procesamiento en
paralelo era plausible y aumentaba el rendimiento en cerca de N veces. Mostraron
problemas al implementar funciones de transferencia pasa banda los cuales
fueron resueltos por un esquema de implementación híbrido.
La paralelización de los procesos es hoy en día uno de los centros de atención
que brindan posibilidades de aumentar el rendimiento. En el trabajo aquí
presentado seguirá una filosofía de paralelización de H(z) para aumentar el
rendimiento del sistema y permitir un mayor control del mismo, como fue mostrado
en [2].
En 1988 una década después de que Intel lanzara su primera DSP al mercado
Walter Demmer y Joachim Speidel de Philips Kommunikations Industrie AG en [3]
propusieron un filtro digital programable por el usuario para aplicaciones de alta
velocidad con el fin de socavar el problema de la sobre carga por procesamiento
de las DSPs de aquel entonces. La propuesta consistía en filtrar la señal por
medio de un conjunto de multiplicadores que realizan su tarea de forma
concurrente. La concepción del sistema de filtrado planteaba una metodología
donde se seleccionaban que multiplicadores se utilizarían y se asignaban cada
uno de sus valores. En tiempo de ejecución, el filtro recibía la señal de entrada,
esta pasaba a través de cada uno de los multiplicadores al mismo tiempo, y
finalmente, se sumaban para hallar la salida resultante. El sistema programable
propuesto tenía una gran flexibilidad y permitía implementar una vasta gama de
filtros, puesto que contaba con la posibilidad de implementar hasta 32 coeficientes
totalmente programables por el usuario.
Las técnicas matemáticas para el tratamiento de la función de transferencia de
un filtro digital con el tiempo se han hecho más completas y robustas, lo cual
permitió implementar aplicaciones cada vez más complejas. Particularmente las
técnicas del procesamiento digital de imágenes exigen una carga por
procesamiento bastante alta, por lo cual siempre se ha buscado optimizar dichas
técnicas y mejorar el hardware en el que se implementan para permitirle al sistema
trabajar en tiempo real. En 2004 Jonne Poikonen y Ari Paasio miembros de la
IEEE propusieron un sistema de filtrado de orden superior basado en un modelo
16
de interconexión de “células” no lineales de procesamiento [4], el cual fue
implementado sobre tecnología digital estándar CMOS de 0.18µm. Dicha
realización operaba en paralelo y contemplaba la funcionalidad de orden
estadístico del filtro, haciendo de esta una muy robusta implementación. Los
resultados obtenidos comprobaron el potencial de esta metodología, ya que se
obtuvo una mayor eficiencia en cuanto al uso del espacio de la implementación y
la potencia consumida, además de permitir que dichos procesos operasen
verdaderamente en paralelo para obtener un muy alto rendimiento.
El tratamiento digital de señales, por sus facilidades de diseño e
implementación hace que generalmente, la realización de un filtro se lleve a cabo
en dicho espacio. Sin embargo la implementación de un filtro analógico presenta
ventajas sobre las realizaciones digitales, principalmente que estos no presentan
latencia y por lo tanto no hay pérdida de datos por ella. El problema surge cuando
se desea implementar un filtro IIR (filtrado analógico) en un espacio discreto
(filtrado digital). Para solucionar esta dificultad se tienen muchas técnicas, las más
comunes mapean el sistema en espacio S y lo transforman a espacio Z mediante
transformaciones bilineales. En 2010 Guo Gaizhi, Zhang Pengju, Yu Zongzuo y
Wang Hailong de la Inner Mongolia Normal University implementaron sobre una
DSP TMS320VC54x de la Texas Instruments® un filtro FIR el cual se basó en la
respuesta en frecuencia deseada de su equivalente análogo [5]. El filtro fue
diseñado en Matlab® a partir de los requisitos sobre la respuesta en magnitud. Los
resultados obtenidos concordaron con la simulación hecha en Matlab por lo cual el
método de diseño y realización quedó validado.
Las DSPs han sido una herramienta valiosa para los investigadores desde su
comercialización por las facilidades de diseño, implementación y ejecución de
algoritmos de procesamiento para diferentes propósitos. Sin embargo existen
otras herramientas que pueden llegar a cumplir el mismo objetivo. En 2010 Zhang
Jian, Pan Ling y Ding Quan pertenecientes a la Universidad de Electricidad de
Potencia e Ingeniería de Automatización de Shanghái mostraron en [6] la
implementación de filtros IIR sobre DSPs y tecnología SOPC la cual se basa en
FPGAs. Los resultados experimentales que obtuvieron mostraron que la
implementación con tecnología SOPC posee ventajas sobre las implementaciones
en DSPs bajo un entorno de desarrollo integrado (específicamente Code
Composer Studio [7]). Las principales ventajas observadas fueron controles más
flexibles sobre el sistema, salida de fase continua, alta precisión y confiabilidad. Es
importante resaltar que las ventajas de implementación sobre tecnología SOPC
son producto de que es una realización basada en compuertas por lo cual
básicamente es la implementación directa del algoritmo, mientras que una DSP
17
funciona a base de instrucciones (siguiendo el concepto de una máquina de Turing
[8]), por lo cual una tarea puede tomar más de un ciclo de trabajo.
En la actualidad los sistemas de filtrado avanzados se basan en algoritmos
estadísticos y adaptativos para hacer robusta la aplicación. Los algoritmos
evolutivos presentan un gran potencial para dominar el futuro del procesamiento
puesto que combinan precisión con adaptabilidad y en general, presentan un
mejor resultado con el tiempo que el obtenido implementando algoritmos
tradicionales para el mismo propósito. En 2011 Yao Li y Xuehua Zhang de la
escuela de física de la Universidad de Beihua en China, implementaron un método
de hardware evolutivo para el filtrado mediante variables de estado y usando
amplificadores de transconductancia [9]. Allí en primera instancia determinaban el
tipo, estructura y funcionamiento del filtro de transconductancia, seguido, obtenían
la función de transferencia que este podría implementar y a partir de las técnicas
de diseño de filtros se obtenían los coeficientes correspondientes. A continuación
el desempeño del filtro era evaluado por el algoritmo evolutivo para hallar los
coeficientes de rendimiento y la optimización del sistema, por último
implementaban la optimización en un dispositivo análogo programable. La
supervisión de los coeficientes de rendimiento era llevada a cabo por un módulo
de control y este por medio de un módulo de migración actualizaba el circuito
análogo programable.
18
2. SEÑALES Y SISTEMAS
2.1 Señales Analógicas
Para el estudio de algún fenómeno físico es necesario identificar las
interacciones que este lleva a cabo con el medio en el que se encuentra presente.
El producto de estas interacciones es una señal continua que transporta la
información del fenómeno mediante alguna magnitud física. Las señales son
tratadas como funciones de variables reales permitiendo el estudio de sus
características frecuenciales y la aplicación del cálculo diferencial sobre ellas.
Para el análisis de las señales analógicas una de las herramientas más
importantes es la transformada de Fourier ya que esta describe frecuencialmente
la señal de interés (incluyendo la fase entre los armónicos). Esta transformada
está definida como sigue:
F𝑓(𝑡) = 𝐹(𝑤) = ∫ 𝑓(𝑡)𝑒−𝑗𝑤𝑡𝑑𝑡 ∞
−∞ (1)
2.2 Señales discretas
Las señales en tiempo discreto son aquellas funciones de variable
independiente entera, donde, los intervalos entre los puntos del eje independiente
de la función no se encuentran definidos [11]. Las señales en tiempo discreto
pueden ser abstraídas como señales continuas muestreadas y, este proceso
comúnmente se lleva a cabo mediante la multiplicación de un tren de impulsos por
la función en cuestión como se muestra a continuación:
xp(t) = x(t)p(t) con p(t) = ∑ ∂
∞
n=−∞
(t − nT) (2)
La variable T es llamada tiempo de muestreo puesto que las muestras sobre la
función x(t) se obtienen cada T segundos. Al analizar las repercusiones
producidas en el espectro de la señal x(t) (en espacio w) por la multiplicación del
tren de impulsos, se puede deducir que:
Xp(jw) =1
T∑ X(j(w − kws))
∞
k=−∞
con ws (3)
=2π
T donde ws es llamada frecuencia de muestreo
19
Lo cual conlleva a que Xp(jw) es una función periódica respecto a w. Su
imagen es la superposición periódica del espectro original de la función x(t)
escalada y desplazada. Ahora, suponiendo que x(t) en su espectro posee una
frecuencia wM tal que ∀ |w| > wM → |X(jw)|2 tiende a 0, de allí para que el
espectro de Xp(jw) no posea superposiciones de las réplicas del espectro original
de x(t) entonces se concluye que ws y wM deben cumplir la siguiente inecuación:
ws > 2wM (4)
La expresión anterior es conocida como el Teorema del muestreo de Nyquist-
Shannon [12] y especifica que para muestrear correctamente una señal, la
frecuencia de muestreo debe ser por lo menos dos veces la frecuencia del máximo
componente espectral de la señal a muestrear.
2.3 Sistemas Físicos
Un sistema físico es un conjunto de elementos, objetos o entidades las cuales
llevan entre sí algún tipo de interacción mediante fenómenos físicos, de tal forma
que dichas interacciones pueden ser modeladas matemáticamente [13]. Los
sistemas de interés en este trabajo son los lineales e invariantes en el tiempo
(LTI), los cuales se pueden describir matemáticamente por su respuesta en
frecuencia y puesto que en general se pueden entender como sistemas de filtrado
frecuencial en adelante se hará referencia a ellos como sistemas de filtrado.
2.4 Sistemas en tiempo continuo
El modelo matemático que describe un sistema de filtrado LTI mediante
ecuaciones diferenciales debe cumplir con la siguiente expresión:
y(n) + a1y(n−1) + a2y
(n−2) + ⋯+ +an−1y′ + any
= b0x(n) + b1x
(n−1) + b2x(n−2) + ⋯+ +bn−1x
′ + bnx (5)
Donde 𝑦 es la salida del sistema y 𝑥 es la entrada al mismo.
Una de las formas que se tiene para describir el comportamiento de sistemas
de filtrado en tiempo continuo es mediante la aplicación de la transformada
unilateral de Laplace [14] sobre la ecuación diferencial que modela el sistema de
interés. Dicha transformada es mostrada a continuación:
20
Lx(t) = X(s) = ∫ x(t)e−stdt∞
0 (6)
Como resultado de la aplicación de la transformada anterior sobre el modelo en
ecuaciones diferenciales de un sistema se llega a una expresión como la
siguiente:
H(s) = Y(s)
X(s)=
A(s)
B(s) (7)
De allí, se tiene que 𝐻(𝑠) es una función racional donde su numerador y
denominador son funciones polinómicas que albergan la información de respuesta
transitoria y estacionaria del sistema respectivamente.
Puesto que los polinomios contenidos por H(s) tienen coeficientes reales es
bien conocido que las posibles raíces que estos posean pueden ser solo de dos
tipos: Reales o complejos conjugados. En la gran parte de la bibliografía referente
a este tema las raíces de los polinomios A(s) y B(s) son llamados ceros y polos
respectivamente, además, también es importante resaltar que la estabilidad de un
sistema LTI viene determinada por la ubicación de sus polos en el plano complejo
puesto que, si todos ellos están ubicados en la parte izquierda del plano, el
sistema se considera estable.
2.5 Sistemas en tiempo discreto
El modelo matemático que describe un sistema de filtrado LTI en tiempo
discreto mediante ecuaciones en diferencia debe cumplir con la siguiente
expresión:
y[n] + a1y[n − 1] + a2y[n − 2] + ⋯+ +an−1y[1] + any[0]
= b0x[n] + b1x[n − 1] + ⋯+ +bn−1x[1] + bnx[0] (8)
Una de las formas que se tiene para describir el comportamiento de sistemas
de filtrado en tiempo discreto es mediante la aplicación de la transformada
unilateral Z sobre la ecuación en diferencias que modela el sistema de interés.
Dicha transformada es mostrada a continuación:
21
Zx[t] = X(z) = ∑ x[t]z−n∞n=0 (9)
Análogamente al aplicar la transformada unilateral Z a la ecuación en
diferencias del sistema se obtiene una expresión de la forma:
H(z) = Y(z)
X(z)=
A(z)
B(z) (10)
H(z) Se define como la función de transferencia del sistema en tiempo discreto
y sus propiedades funcionales son análogas a la función de transferencia H(s)
que se mostró anteriormente.
2.6 Relación entre sistemas en tiempo continuo y discreto
Los sistemas en tiempo continuo se pueden llevar a espacio discreto para
realizar allí su implementación. Hay diversas formas de hacerlo, el común es
emplear una transformación de modelo y cambiar el espacio donde se encuentra.
H(s) → Tranformacion → H(z)
Estas transformaciones comúnmente intentan hallar una relación que conserve
las propiedades de H(s) pero en el espacio resultante. La Tabla 1 muestra algunas
de las transformaciones más utilizadas.
Tabla 1. Transformaciones para llevar sistemas en tiempo continúo a tiempo discreto.
# Nombre Formula
1
Transformación Bilineal. s ←
2
T
z − 1
z + 1
2
Retención de orden cero
u(t) = u[k],
kTz ≤ t ≤ (k + 1)Tz
22
3
Polos y ceros equivalentes
zi = esiTs
2.7 Bancos de filtros [15]
Es un conjunto de filtros organizados en paralelo que comparten una entrada
común y poseen una única salida consolidada. En general se emplean para el
procesamiento específico de las diferentes bandas de interés en la señal que se
desee estudiar.
Un banco de filtros consta de tres fases, fase de Análisis, fase de
Procesamiento y fase de Síntesis. A continuación se describe el objetivo de cada
fase:
o Fase de Análisis:
En esta fase se lleva a cabo la separación espectral de las diferentes
componentes de interés de la señal a estudiar.
o Fase de Procesamiento:
Esta fase tiene como propósito darle el espacio al implementador del
proyecto para llevar a cabo cualquier tipo de procesamiento que
desee con las muestras producto de la fase de análisis.
o Fase de síntesis:
El objetivo de esta fase es reconstruir a partir de los resultados de la
fase de procesamiento la señal de salida que se tenga como
objetivo.
23
3. PROCESADORES DIGITALES DE SEÑALES Y HERRAMIENTAS
ARQUITECTÓNICAS
3.1 Procesadores digitales de señales (DSPs)
Las DSPs son dispositivos cuyos componentes (tanto hardware como
software) han sido optimizados para realizar tareas numéricas complejas donde se
requiera un alto desempeño y eficiencia. Generalmente por sus características son
usadas en aplicaciones de alto desempeño para el procesamiento en tiempo real
de señales analógicas como lo son: el reconocimiento de voz, telefonía celular,
procesamiento de audio y video, etc.
3.2 ADSP Blackfin-533 [16]
Es una DSP fabricada por Analog Devices® y es sobre la cual el presente
trabajo se implementó. La ADSP-BF533 es una DSP de propósito general la cual
pertenece a la familia de productos Blackfin® de Analog Devices y cuenta con una
arquitectura de micro señal (MSA), amplia gama de voltajes de funcionamiento,
regulador de voltaje programable, diversos módulos que facilitan las aplicaciones
de multimedia, entre muchos otros componentes.
La ADSP-BF533 es un dispositivo para aplicaciones que requieran de un alto
performance ya que cuenta con una arquitectura de cache flexible, sistema de
acceso directo a memoria (DMA), velocidad de reloj principal de hasta 600MHz,
módulo para el manejo dinámico del consumo de energía (DPM) produciendo así
un dispositivo flexible, óptimo para una gran gama de aplicaciones.
Entre los componentes que permiten a la ADSP BF-533 ser un dispositivo de
alto rendimiento están: Dos multiplicadores acumulativos (MAC), arquitectura de
procesamiento de múltiples datos (SIMD), componentes embebidos para el
procesamiento multimedia, conjunto de instrucciones reducido (RISC) entre otras
características. Los periféricos principales de dicho dispositivo posee son: Timmer
de vigilancia (WatchDog), interfaz para periférico serial (SPI), sistema para
periféricos paralelos (PPI), transmisor/receptor asíncrono universal (UART), reloj
de tiempo real (RTC), puertos entrada/salida de propósito general además de un
controlador de acceso directo a memoria (DMA) cuya función es transmitir la
información de dispositivos externos a las memorias internas sin intervención del
procesador.
24
La ADSP BF-533 también cuenta con un códec para el procesamiento de audio
cuya referencia es AD1836 y el cual es el usado para la implementación de filtros
mediante la metodología presentada en este trabajo. Las principales
características del códec AD1836 se mencionan a continuación:
- Sistema de audio de 5v multicanal.
- Soporte para procesar las señales a 24 bits con frecuencia de muestreo d
e 96 KHz
- Salida diferencial para un óptimo desempeño.
- ADCs: -92 dB THD + N, 105 dB SNR con rango dinámico.
- DACs: -95 dB THD + N, 108 dB SNR con rango dinámico.
- Ganancia programable para las entradas ADC.
- Puerto de datos serial flexible con justificación a la izquierda o derecha,
compatible con I2S y modo para puerto serial de DSP.
- Frecuencias de muestreo compatibles: 32 KHz, 44 KHz, 48 KHz y 96 KHz.
Es necesario mencionar que se pueden manejar tres salidas y dos entradas
simultáneamente a una frecuencia de muestreo máxima de 48 KHz, ya que para
trabajar a 96 KHz solo se pueden manejar como máximo, dos entradas y dos
salidas simultáneamente. El hardware de la ADSP BF-533 es mostrado en la
Figura 1 donde además, se puede observar el códec, las entradas de audio a la
DSP y las respectivas salidas.
Figura 1. Hardware de la ADSP-BF533.
25
3.3 Arquitectura de BUS y Memoria de la ADSP BF-533
Figura 2. Modelo de memoria Jerárquico de la ADSP BF-533.
La arquitectura en la que se basa la ADSP-BF533 es la Harvard con algunas
modificaciones sobre el manejo de la memoria que la hacen más eficiente. La
memoria de la BF533 es de uso compartido entre la memoria de datos y la de
programa donde se sigue el patrón jerárquico mostrado en la Figura 2 cuyo orden
es de más lenta (Memoria externa L2) a la más rápida (memoria interna L1). La
arquitectura general de la memoria de la ADSP-BF533 es mostrada en la Figura 3
allí se puede observar que el procesador cuenta con 4 buses principales:
- Un bus de 64 bits para las instrucciones.
- Dos buses de 32 bits para el cargue simultaneo de datos.
- Un bus de 32 bits para el almacenamiento de datos.
Además, es importante hacer hincapié en que el reloj de funcionamiento de la
memoria L1 es el mismo reloj del procesador principal (CCLK).
Figura 3. Arquitectura de memoria de la ADSP-BF533.
26
Figura 4. Mapa de la memoria de la ADSP-BF5xx.
El mapa de memoria mostrado en la Figura 4 permite observar que la DSP
tiene:
- 80 Kb de memoria de programa 0xFFA0 0000 hasta 0xFFA1 3FFF
- 32 Kb de memoria de datos (Banco A) 0xFF80 0000 hasta 0xFF80 8000
- 32 Kb de memoria de datos (Banco B) 0xFF90 0000 hasta 0xFF90 8000
- 4 K bytes de memoria de block de notas para almacenar datos generales,
como la pila.
Con lo cual claramente se tiene un total de 148 KBytes de memoria interna
disponibles. La memoria interna es calificada con L1 y algunas secciones de ella
pueden ser configuradas como memoria estática de acceso aleatorio (SRAM).
La BF533 también posee una memoria on-chip (ROM) de arranque para
inicializar la tarjeta y permitir cargar memoria de programas que se encuentren en
memorias externas a la L1.
27
3.4 Archivos de registros y unidades de procesamiento de Hardware
En la Figura 5 se puede observar los tres principales componentes del núcleo
del procesador ADSP-BF533 los cuales son: Unidad aritmética de dirección,
unidad de control y la unidad aritmética de datos.
Unidad aritmética de dirección
Componentes:
- Dos generadores de dirección de datos (DAG1 Y DAG2).
- Seis registros de puntero de dirección de propósito general de 32 bits.
- Un puntero de marca de 32 bits (FP).
- Un puntero de pila de 32 bits (SP).
- Un set de registros generadores de dirección de datos de 32 bits.
La función principal de la unidad aritmética de dirección consiste en generar
las direcciones de acceso a la memoria del procesador mediante sus 6
punteros de direcciones.
Figura 5. Núcleo de la ADSP-BF533 [17].
28
Unidad de control
Componentes:
- Secuenciador de las instrucciones de programa para controlar la ejecución
de instrucciones donde se incluyen también las instrucciones de alineación y
de decodificación.
- Buffer de ciclo para dar soporte a los bucles de cero sobrecarga.
Unidad aritmética de datos
Componentes:
- Dos multiplicadores de 16 bits.
- Dos acumuladores de 40 bits (ACC0 Y ACC1).
- Dos unidades aritmético lógicas ALUs de 40 bits.
- Cuatro unidades aritmético lógicas de video de 8.
- Un Shifter de 40 bits.
- Ocho registros de datos de 32 bits (R0 A R7).
Estos componentes son los más utilizados por las aplicaciones industriales
puesto que es aquí donde se lleva a cabo el procesamiento de los datos en sus
diferentes formatos. La información para el cálculo llega desde los buses de
datos hacia la unidad aritmética de datos y cuando estos son procesados los
resultados se trasladan a los registros de datos antes de transferirlos a la
memoria a través de los buses de datos.
3.5 Periféricos Básicos.
En la ADSP-BF533 se tiene bastante flexibilidad en cuanto al manejo de
periféricos puesto que tiene a su disposición bastantes herramientas para el
manejo de ellos. En la Figura 3 se puede observar algunos de los periféricos que
aquí son listados:
- Un temporizador watchdog controlado por el reloj del sistema (SCLK).
- Un reloj en tiempo real el cual funciona como un reloj digital para el
procesador. Puede funcionar como alarma para la activación de la
tarjeta.
- Tres temporizadores de propósito general son configurados como PWM,
captura de ancho y período, y contador de eventos externos.
- 16 flags bidireccionales programables de propósito general (PF0-PF15).
29
- Un receptor/transmisor universal asíncrono (UART), tiene una velocidad
de transmisión máxima de hasta SCLK/16.
- Dos puertos seriales asíncronos (SPORT0 Y SPORT1) que permiten
comunicación serial de alta velocidad, con un velocidad máxima
SCLK/2.
- Un SPI para comunicación serial de alta velocidad, de máximo SCLK/4.
- Un PPI permite el manejo de un bus paralelo programable con una
velocidad máxima de SCLK/2. Se suele usar como conversor de datos
de alta velocidad.
- EBIU para permitir interfaz glueless1 con memorias externas. Tres de los
buses internos de 16 bits están unidos a la EBIU.
- El bus de acceso externo (EAB) está controlado por la memoria de
núcleo para acceder a la memoria externa.
- El bus de acceso periférico (PAB) se utiliza cuando se necesita ingresar
a los registros de la EBIU.
3.6 Entorno de Desarrollo Integrado. La programación de dispositivos digitales por lo general se lleva a cabo sobre
un ambiente de desarrollo que permite implementar, mediante codificación,
algoritmos que realicen las tareas que se necesiten. Estos ambientes son
llamados IDEs (de sus siglas en inglés Integrated Development Environment) y
cumplen un rol fundamental al momento de programar todo tipo de dispositivos
puesto que ellos contienen todas las herramientas necesarias para la codificación
y la puesta a punto de la realización computacional de los algoritmos a
implementar. Una de las características más importantes de un IDE es su
depurador, puesto que este es el que le permite ver al desarrollador la ejecución
paso a paso del algoritmo que se implementó para la corrección de errores y
optimización de código en el algoritmo.
3.7 Visual DSP ++ 5.0 [18] El IDE utilizado para la implementación de los algoritmos en este trabajo fue el
Visual DSP++ 5.0® el cual es un entorno para el desarrollo de proyectos y
1 Término usado cuando no es necesario tener un circuito entre dos módulos ICs o bloques de hardware.
30
aplicaciones de cualquier tipo cuyas implementaciones se realizan sobre DSPs.
Visual DSP++ permite el manejo de los siguientes lenguajes de programación:
- C
- C++
- ASM (Assembly)
Además, Visual DSP++ tiene un conjunto de herramientas que vienen
determinadas por el tipo de licencia que se libere. Entre las básicas se encuentran
el compilador C/C++, un simulador de nivel de ciclo, un archivador, un depurador,
una librería para tiempo de ejecución y ciertas librerías para un desarrollo ágil.
Visual DSP++ contiene algunas características que son importantes resaltar
puesto que son de vital interés para el manejo a bajo nivel de las tarjetas. El
VisualDSP++ Kernel (VDK) es una herramienta que permite la gestión de los
recursos propios de las tarjetas, con lo cual las implementaciones que requieran
un alto nivel de eficiencia puedan ser implementadas minuciosamente.
Figura 6. Visual DSP++ 5.0.
3.8 Manejo de proyectos sobre Visual DSP ++ 5.0
Para realizar depuración sobre el código que se ha escrito se deben realizar los
siguientes pasos:
- Establecer conexión con la “Session” de depuración.
- Creación del proyecto y adición del código al mismo.
- Ejecutar en modo Debug (depuración).
31
El primer paso es establecer conexión con la tarjeta en la cual se correrá el
programa, para ello se selecciona la pestaña de SESSION y luego a NEW
SESSION, al darle clic allí aparecerá una ventana como la que se muestra en la
Figura 7 donde se debe diligenciar los siguientes aspectos de la Sesión:
- DSP que se desea utilizar para la implementación del código.
- La conexión a utilizar, conexión física a la tarjeta o un emulador de dicha
tarjeta.
- Plataforma para la depuración del código.
Figura 7. Ventana para la creación/edición de Sesiones en Visual DSP++.
Luego de crear efectivamente la sesión se debe crear un proyecto y para esto
se debe de ir a la pestaña FILE y a la opción de NEW/PROJECT donde aparecerá
una venta para la configuración del nuevo proyecto como se muestra en la Figura
8, allí se deben configurar los siguientes elementos para la creación del proyecto:
Tipos de proyecto:
- Standard Aplication: Son los proyectos estándar que no requieren
herramientas avanzadas del IDE.
- Library: Son proyectos cuyo resultado son librerías que se pueden
utilizar para otros proyectos.
- LwIP Ethernet A: Son proyectos cuyo fin son aplicaciones de red.
32
- VDK Aplication: Son proyectos en cuya arquitectura la librería para el
manejo del Kernel VDK es fundamental.
Procesador a utilizar: Se debe seleccionar el procesador ó tarjeta sobre la
cual se va a implementar el proyecto.
Configuraciones de la Aplicación: Dependiendo del tipo de proyecto las
configuraciones varían.
Añadir Código inicial: Al momento de la creación del proyecto Visual DSP++
cuenta con la opción de registrar un código inicial al proyecto para el
entendimiento del funcionamiento del tipo de proyecto.
Figura 8. Cuadro para la configuración de nuevo proyecto en Visual DSP++.
3.9 Depuración sobre Visual DSP ++ 5.0 Para empezar una sesión de depuración sobre el proyecto se debe seleccionar
la opción de Debug y el proyecto a debuguear acto seguido se debe dar clic en
REBUILT ALL como se muestra en la Figura 9. Al darle clic en REBUILT ALL se le
indica a Visual DSP++ que compile los archivos del proyecto y los ingrese en la
memoria de programa del procesador seleccionado para así permitirle al usuario la
ejecución de la aplicación. Las instrucciones que permite el Visual DSP++ para
realizar la depuración son las mostradas en la Tabla 2.
33
Tabla 2. Herramientas para la depuración de código en un proyecto sobre Visual DSP++.
Comando Descripción
Run
Corre un programa ejecutable. El programa corre hasta que un evento lo detenga, tal como un
punto de corte o intervención del usuario. Cuando la ejecución del programa se detiene, todas las
ventanas se actualizan para la actual dirección y valor
Halt
Detiene la ejecución del programa. Todas las ventanas se actualizan después que el procesador se
detiene. Los valores que cambiaron se destacan, y el estado de la barra despliega la dirección
donde el programa se detiene.
Run to Cursor Corre el programa a la línea donde está el cursor. Se puede colocar el cursor en la ventana de
edición y en la ventana de Disassembly.
Step Over
(Código C/C++ Solo en la ventana de edición) Pasos simples hacia adelante través de las
instrucciones de programa. Si la línea de origen llama a una función, la función es ejecutada
completamente, sin efectuar paso a paso las instrucciones de programa.
Step Into (Ventana de edición o ventana de Diassembly) Pasos simples a través del programa una instrucción
C/C++ o assembly a la vez. Las funciones encontradas son ingresadas.
Step Out Of
(Código C/C++ solo en la ventana de edición) Realiza pasos múltiples hasta que la función actual
retorna su llamado, y se detiene en la instrucción que sigue inmediatamente el llamado de la
función.
Figura 9. Configuración para la depuración de un proyecto.
34
Para evaluar el comportamiento de los datos que ingresan a la DSP, Visual
DSP++ cuenta con diversas herramientas, entre las más útiles se encuentra el
visor gráfico para el monitoreo del comportamiento de variables. Para configurar
esta herramienta se puede realizar el proceso mostrado en la Figura 10 y para su
uso se debe configurar el lugar preciso en la memoria donde se encuentran las
variables a utilizar para el gráfico, indicar el tipo de variable, el número de
muestras entre otros datos que se solicitan.
Figura 10. Procedimiento para la utilización de la herramienta de monitoreo gráfico de variables.
3.10 Herramientas de medición: Tarjetas de Adquisición de datos. Las tarjetas de adquisición de datos son dispositivos diseñados para la
adquisición de medidas sobre los fenómenos físicos contando entre sus
principales características la fidelidad, resolución y fiabilidad en las medidas
realizadas. Generalmente estas tarjetas son denominadas DAQs (Data Acquisition
35
Card) y son los instrumentos usados para convertir las medidas analógicas en
digitales y poder llevar dichas medidas a sistemas de procesamiento.
3.11 NI MyDAQ [19]. La tarjeta de adquisición NI MyDAQ es un dispositivo diseñado y distribuido por
National Instruments®. Esta es una DAQ de bajo costo, portable y cuyo software
instrumental esta implementado sobre LabVIEW. Este dispositivo permite la
medición de señales del mundo real y es ideal para explorar la electrónica de la
toma de medidas provenientes de sensores. La tarjeta físicamente se muestra en
la Figura 11.
Figura 11. Tarjeta NI MyDAQ.
3.12 Estructura física de la NI MyDAQ. La tarjeta NI MyDAQ cuenta con un lateral de 20 pines donde se encuentran
ubicadas las diferentes conexiones que esta tarjeta posee, además cuenta con 2
conectores de audio mono para la entrada y salida de datos analógicos (Figura
12). La NI MyDAQ cuenta también, con tres conectores para el uso del multímetro
digital y estos se pueden observar en la Figura 13.
Figura 12. Lateral de la NI MyDAQ.
36
Figura 13. Conectores de la NI MyDAQ para el uso del multímetro digital. (1) Medición de: Resistencia, voltaje, diodo y continuidad. (2) Medición de corriente.
3.13 Componentes del Hardware de la NI MyDAQ. La NI MyDAQ cuenta con los componentes mostrados en la Figura 14, donde
se deben resaltar los siguientes:
- 8 Canales de Entrada/Salida Digitales.
- 8 Canales de Entrada Análogas.
- 2 Canales de Salida Análogas.
- Multímetro Digital.
- Reguladores de Tensión.
Figura 14. Diagrama de bloques de los Componentes de la NI MyDAQ.
37
Canales de Entrada/Salida Digitales (Digital I/O)
La NI MyDAQ cuenta con 8 canales de propósito general para la toma o
el envió de datos en forma de pulsos digitales. Estos canales están
diseñados para leer a 3.3v y 5v y es importante resaltar que la salida digital
de la NI MyDAQ no es compatible con 5v CMOS.
Canales de Entrada Análogas (AI)
La NI MyDAQ cuenta con dos entradas analógicas las cuales pueden
ser configuradas como canales de propósito general de alta impedancia.
Las entradas analógicas son multiplexadas y llevadas a un único ADC el
cual cuenta con un rango de funcionamiento de ±10V. Las medidas
extraídas del ADC son entregadas hasta 200 kS/s por canal y por esto la NI
MyDAQ permite tener aplicaciones como Osciloscopio, Analizador dinámico
de señales y Analizador de Bode.
Canales de Salida Análogas (AO)
La NI MyDAQ cuenta con dos salidas analógicas, estas pueden ser
configuradas como salidas de propósito general o salidas de audio. Ambos
canales generan señales en el rango de ±10V. Si estos son configurados
como salidas de audio se comportan como canal izquierdo y derecho
respectivamente. Las salidas en estos canales son actualizadas hasta 200
kS/s por cada canal lo cual permite que este dispositivo sea un generador
de señales y analizador de bode adecuado para la instrumentación a bajas
frecuencias.
Multímetro Digital
La NI MyDAQ provee funcionalidades para la medición de voltaje (AC y
DC), corriente (DC y AC), resistencia, y voltaje de caída de diodo. Este
multímetro cuenta con un software el cual está en el CD de instalación y
permite realizar las diferentes medidas mencionadas mediante una GUI en
la computadora.
Reguladores de Tensión
La NI MyDAQ permite obtener fuentes de voltaje constante entre ±15V
las cuales pueden ser usadas como elemento activo para componentes
analógicos como amplificadores operacionales y reguladores lineales. La NI
MyDAQ cuenta también con otra fuente de 5V la cual suele utilizarse como
38
alimentación para dispositivos lógicos. La potencia total que puede soportar
estas fuentes es de 500 mW donde el cálculo de potencia se lleva a cabo
mediante el teorema de Joule.
3.14 Funcionalidades pre-programadas para la NI MyDAQ. La NI MyDAQ cuenta con las siguientes funcionalidades pre-programadas para
el uso inmediato luego de su instalación:
Multímetro Digital (DMM)
El NI ELVISmx multímetro digital permite la medición de variables
eléctricas mediante el uso del harware de la NI MyDAQ y el procesamiento
mediante LabVIEW de la GUI NI ELVISmx. Los parámetros de medición del
multímetro digital son mostrados a continuación:
- Rangos de voltaje DC: 60V, 20V, 2V y 200mV.
- Rangos de voltaje AC: 20V, 2V y 200 mV.
- Rangos de corriente DC: 1A, 200mA y 20mA.
- Rangos de corriente AC: 1A, 200mA y 20mA.
- Rangos de Resistencia: 20MΩ, 2MΩ, 200kΩ, 20kΩ, 2kΩ y 200Ω.
- Rango de diodo: 2v.
- Resolución (número de cifras significativas para mostrar): 3.5.
Figura 15. Interfaz para el uso del multímetro digital de la NI MyDAQ.
39
Osciloscopio digital
El osciloscopio NI ELVISmx es un medidor de voltaje en tiempo real para el
análisis gráfico de estos datos. En Figura 16 se puede observar la GUI para
esta funcionalidad. Cuenta con los siguientes parámetros de medición:
- Canal Origen: canal por donde se medirá la señal de interés.
Puede ser alguno de los siguientes canales AI0, AI1,
AudioInputLeft o AudioInputRigth.
- Escala de Voltaje: para canales AI- 5V, 2V, 1V, 500mV, 200mV,
100mV, 20mV y para entradas de Audio- 1V, 500mV, 200mV,
100mV, 50mV, 20mV.
- Tiempo de base: Es permitido configurar valores entre el rango de
200ms hasta 5µs.
- Trigger: Se permite el inmediato y el edgetrigger.
Figura 16. Interfaz para el uso de la funcionalidad: Osciloscopio digital de la NI MyDAQ.
Generador de Funciones (FGEN)
El generador de señales de NI ELVISmx genera funciones de onda
estándar donde sus parámetros pueden ser variados en tiempo de
ejecución. Puede generar ondas sinusoidales, cuadradas o triangulares.
Los parámetros de esta funcionalidad pre programada son:
40
- Canal de salida: El canal AO0 o AO1.
- Rango de frecuencia: 0.2Hz hasta 20kHz.
Figura 17. Interfaz para el uso de la funcionalidad: Generador de funciones de la NI MyDAQ.
Analizador de Bode (BODE)
El analizador de Bode de la NI ELVISmx es quizás la herramienta más útil
que la NI MyDAQ tiene pre-programada puesto que esta funcionalidad
permite ver rápidamente la respuesta en frecuencia para la magnitud y la
fase de un sistema. La interfaz de esta funcionalidad es mostrada en la
Figura 18. Los parámetros de medición son:
- La medida de la entrada al sistema se realiza sobre el canal AI 0.
- La medida de la salida al sistema se realiza sobre el canal AI 1.
- La señal de estímulo se produce en el canal AO 0.
- El rango de frecuencias de análisis es de 1Hz hasta 20kHz.
41
Figura 18. Interfaz para el uso de la funcionalidad: Analizador de Bode de la NI MyDAQ.
Analizador dinámico de señales (DSA)
El analizador dinámico de la NI ELVISmx realiza una transformación al
dominio de la frecuencia de la forma de onda medida en el canal AI o en el
de audio. El DSA también permite aplicar distintas ventanas para el
muestreo de la señal. La GUI mostrada en la Figura 19 es la interfaz para el
manejo del DSA de la NI ELVISmx.
Figura 19. Interfaz para el uso de la funcionalidad: Analizador dinámico de señales de la NI MyDAQ.
42
4. METODOLOGÍA
En la actualidad la paralelización de los sistemas de procesamiento es una
herramienta de implementación fundamental cuando se requiere de
procedimientos óptimos para el procesamiento de la información. De allí, la
metodología de implementación aquí mostrada pretende brindar un camino fácil y
eficiente para la implementación de filtros LTI digitales y en paralelo.
4.1 Perspectiva General y Consideraciones Iniciales La metodología plantea la subdivisión del sistema a implementar en
subsistemas los cuales puedan reconstruir efectivamente el sistema original y
además puedan ser implementados en paralelo de forma sencilla por el
desarrollador del proyecto.
La reconstrucción del sistema puede realizarse por medio de una sumatoria,
productoria o cualquier tipo de operación que reconstruya efectivamente el
sistema original. En este proyecto la implementación de la metodología se realizó
a base de sumatorias.
El desarrollo metodológico mostrado a continuación parte de la función de
transferencia del sistema, sin embargo, el concepto general de reconstrucción
mostrado, puede ser aplicado a cualquier método de diseño que reduzca el
sistema original a subsistemas de orden uno (polos reales) o dos (polos complejos
conjugados).
4.2 Base Teórica
Para la implementación de un sistema LTI digital mediante su reconstrucción a
partir de subsistemas se puede partir de su función de transferencia como sigue:
H(z) =A(z)
B(z) (11)
De allí reescribiendo se puede reescribir H(z) como:
H(z) =A0(z)
B0(z)+
A1(z)
B1(z)+ ⋯+
An(𝑧)
Bn(z)=
A(z)
B(z)
H(𝑧) = 𝐺0(𝑧) + 𝐺1(𝑧) + ⋯+ 𝐺𝑛(𝑧)
H(𝑧) = ∑𝐺𝑖(𝑧)
𝑛
𝑖=0
(12)
43
Generalizando se puede escribir la expresión anterior como:
H(𝑧) = ∑𝐻i(z)
n
i=0
(13)
Con
𝐻i(z) = CiGi(z) (14)
Los coeficientes de cada subsistema Ci son llamados coeficientes de flexibilidad
puesto que permiten parametrizar la implementación del filtro y moldear la
respuesta en frecuencia de la reconstrucción en función de la ponderación de
cada subsistema. A partir del análisis de las expresiones anteriores se hace
importante resaltar algunas observaciones:
- H(z) representa al sistema en el espacio Z y estará dada por la
sumatoria (13) donde, se puede considerar cada 𝐻i como un
subsistema y a H(z) como una superposición de las respuestas
en frecuencia de estos subsistemas.
- Ya que H(z) es un sistema LTI propio y teniendo en mente que
todos los coeficientes de los polinomios encontrados en H(z) son
reales, entonces es claro que cada 𝐻i, al reducir lo máximo
posible H(z), debe ser un subsistema de primer o segundo orden,
donde este último tendrá polos conjugados complejos.
- El buffer que contenga los datos de la muestra actual de entrada
y las anteriores dos es común a todos los subsistemas 𝐻i.
- Cada subsistema 𝐻i deberá de tener un buffer de salida único
para la muestra actual y las dos anteriores de dicho subsistema
puesto que estos son independientes entre sí.
- La muestra actual de salida del filtro está determinada por la
combinación lineal de todas las salidas de los subsistemas y
además no se tendrá la obligación de almacenar estas muestras
en un buffer.
4.3 Definiciones y Convenciones
A continuación se describen algunas definiciones y convenciones que se usarán
para describir los algoritmos de funcionamiento más adelante:
44
- Número de Subsistemas (M): Es el número de subsistemas en los
cuales el sistema original se ha dividido.
- Matriz de Buferes de salida (MBS): Matriz de tamaño M x 2 donde
se almacenan las muestras de salida de los diferentes
subsistemas.
- Buffer de entrada (BE): Arreglo de 3 posiciones el cual contiene la
muestra actual y las dos anteriores de entrada.
- Matriz de Coeficientes (MCO): Es la matriz donde se almacenan
todos los coeficientes de los subsistemas. Es de tamaño 2M x 3.
- Conjunto de Órdenes (CO): Es un arreglo de M posiciones que
contiene el orden de cada uno de los sistemas contenidos en la
matriz de coeficientes.
- Constante espectral (CE): Es la constante que determina la
ganancia independiente a la frecuencia obtenida por el proceso
de subdivisión.
4.4 Algoritmos y Flujos En la Figura 20 se muestra el proceso general que el proyecto debe de realizar
para implementar cualquier sistema LTI mediante la metodología propuesta, cabe
aclarar que el proceso de subdivisión del sistema a implementar es realizado por
medio del algoritmo de fracciones parciales.
Figura 20. Algoritmo general de la implementación sobre la DSP.
45
A continuación se describirá cada uno de estos pasos enfocados a la
conceptualización teórica. Luego se verán los mismos con un enfoque hacia la
codificación realizada sobre el Visual DSP++.
1. Inicialización de Variables Fundamentales
Se deben inicializar las variables donde se almacenarán los diferentes
valores que son fundamentales para el funcionamiento de la metodología.
A continuación se listan algunas de estas variables:
- 4 Variables para el manejo de los 4 canales de entrada.
- 4 Variables para el manejo de los 4 canales de salida.
- Registro de configuración del códec 1836.
- Buffer de tres posiciones para las muestras de entrada.
- Matriz de Buferes de M x 2 posiciones para las muestras de
salida de cada subsistema.
- Matriz de coeficientes de 2M x 3 para almacenar los coeficientes
de los subsistemas resultantes.
- Polinomio de numerador del sistema original a implementar.
- Polinomio de denominador del sistema original a implementar.
2. Lectura del sistema a Implementar:
La metodología debe poder conocer los coeficientes de la función de
transferencia del filtro original para cumplir con su objetivo de
reconstruirlo.
3. Hallar los polos del sistema a implementar:
Para hallar las raíces del denominador del sistema se puede implementar
cualquier método numérico que pueda hallar raíces reales y complejas
sobre un polinomio de coeficientes reales.
En general se tendrá un polinomio 𝐵(𝑧) de orden N por lo cual se puede
inferir que habrán:
46
𝑁1 𝑃𝑅1, 𝑃𝑅2, … , 𝑃𝑅𝑁1 𝑝𝑜𝑙𝑜𝑠 𝑟𝑒𝑎𝑙𝑒𝑠, 𝑦
𝑁2 𝑃𝐶1, 𝑃𝐶2, … , 𝑃𝐶𝑁2 𝑝𝑜𝑙𝑜𝑠 𝑐𝑜𝑚𝑝𝑙𝑒𝑗𝑜𝑠 𝑐𝑜𝑛𝑗𝑢𝑔𝑎𝑑𝑜𝑠
𝐶𝑜𝑛 𝑁 = 𝑁1 + 2𝑁2
Por lo cual 𝐵(𝑧) se puede ver como:
𝐵(𝑧) = (𝑧 − 𝑃𝑅1)(𝑧 − 𝑃𝑅2)… (𝑧 − 𝑃𝑅𝑁1)(𝑧 − 𝑃𝐶1) ∗ (𝑧 − 𝑃𝐶1
)(𝑧 − 𝑃𝐶2)(𝑧 − 𝑃𝐶2 ) … (𝑧 − 𝑃𝐶𝑁2)(𝑧 − 𝑃𝐶𝑁2
) (15)
4. Plantear y resolver el sistema de ecuaciones:
El objetivo de este paso es partir de la función de transferencia original y
dar como resultado los coeficientes de cada uno de los subsistemas, esto
mediante la aplicación de fracciones parciales al sistema original.
Partiendo de:
𝐻(𝑧) = 𝐴(𝑧)𝐵(𝑧)⁄
Entonces se puede observar que la reconstrucción de 𝐻(𝑧) viene dada
por:
𝐻(𝑧) = 𝐻𝑅(𝑧) + 𝐻𝐶(𝑧) + 𝐶𝐸 (16)
Con:
𝐻𝑅(𝑧) = 𝐶𝑅10𝑧 + 𝐶𝑅20
𝑧 − 𝑃𝑅1+
𝐶𝑅11𝑧 + 𝐶𝑅21
𝑧 − 𝑃𝑅2+ ⋯+
𝐶𝑅1𝑁1𝑧 + 𝐶𝑅2𝑁1
𝑧 − 𝑃𝑅𝑁1 (17),
𝐻𝐶(𝑧) = 𝐶𝐶00𝑧
2 + 𝐶𝐶10𝑧 + 𝐶𝐶20
(𝑧 − 𝑃𝐶1)(𝑧 − 𝑃𝐶1 )
+𝐶𝐶01𝑧
2 + 𝐶𝐶11𝑧 + 𝐶𝐶21
(𝑧 − 𝑃𝐶2)(𝑧 − 𝑃𝐶2 )
+ ⋯+𝐶𝐶0𝑁2𝑧
2 + 𝐶𝐶1𝑁2𝑧 + 𝐶𝐶2𝑁2
(𝑧 − 𝑃𝐶𝑁2)(𝑧 − 𝑃𝐶𝑁2 )
(18)
Y CE es la constante espectral, cuyo valor es el coeficiente de zN en el
numerador del sistema original.
Se debe observar que las constantes 𝐶𝑅𝑖𝐽 y 𝐶𝐶𝑖𝐽 son las incógnitas a
determinar para poder implementar cada uno de los subsistemas puesto
que los coeficientes de los denominadores 𝑃𝑅𝑖 y 𝑃𝐶𝑖 son las raíces del
denominador de la función de transferencia.
El planteamiento del sistema de ecuaciones se realiza mediante el
siguiente concepto de igualdad de polinomios:
47
𝑃𝑎𝑟𝑎 𝑞𝑢𝑒 𝑑𝑜𝑠 𝑝𝑜𝑙𝑖𝑛𝑜𝑚𝑖𝑜𝑠 𝑃(𝑥) 𝑦 𝑄(𝑥) 𝑑𝑒 ó𝑟𝑑𝑒𝑛𝑒𝑠 𝑁 𝑦 𝑀 𝑟𝑒𝑠𝑝𝑒𝑐𝑡𝑖𝑣𝑎𝑚𝑒𝑛𝑡𝑒
𝑠𝑒𝑎𝑛 𝑖𝑔𝑢𝑎𝑙𝑒𝑠 𝑒𝑛𝑡𝑜𝑛𝑐𝑒𝑠 𝑠𝑒 𝑑𝑒𝑏𝑒 𝑐𝑢𝑚𝑝𝑙𝑖𝑟 𝑞𝑢𝑒: 𝑃0 = 𝑄0, 𝑃1 = 𝑄1, … , 𝑃𝑁 = 𝑄𝑀, 𝑑𝑜𝑛𝑑𝑒 𝑁 = 𝑀,𝑄𝑖 𝑦 𝑃𝑖 𝑠𝑜𝑛 𝑙𝑜𝑠 𝑐𝑜𝑒𝑓𝑖𝑐𝑖𝑒𝑛𝑡𝑒𝑠
𝑑𝑒𝑙 𝑡𝑒𝑟𝑚𝑖𝑛𝑜 𝑑𝑒 𝑋𝑁−𝑖 𝑦 𝑋𝑀−𝑖 𝑟𝑒𝑠𝑝𝑒𝑐𝑡𝑖𝑣𝑎𝑚𝑒𝑛𝑡𝑒.
Partiendo de la ecuación (16) y teniendo en cuenta que el denominador
de la parte izquierda y derecha son iguales (por el procedimiento que se
siguió) entonces el sistema a resolver es producto de la igualdad en los
numeradores de dicha expresión.
El planteamiento que se realiza para el sistema es el siguiente:
𝑀 =
[
𝑀00 𝑀01
𝑀10 𝑀11⋯
𝑀0,𝑁−1 𝑀0𝑁
𝑀1,𝑁−1 𝑀1,𝑁
⋮ ⋱ ⋮𝑀𝑁−1,0 𝑀𝑁−1,1
𝑀𝑁,0 𝑀𝑁,1⋯
𝑀𝑁−1,𝑁1−1 𝑀𝑁−1,𝑁
𝑀𝑁,𝑁−1 𝑀𝑁,𝑁 ]
= 𝐴 = [
0𝐴1
⋮𝐴𝑁
] (19)
Donde es el vector que contiene todas las incógnitas, 𝐴 el vector con
los coeficientes del numerador del sistema original y 𝑀 la matriz de
relaciones entre las incógnitas y los coeficientes del numerador de la
función de transferencia original. Cabe destacar que siempre se tendrá
para el sistema de ecuaciones 𝐴0 = 0 puesto que este se extrae como la
constante espectral.
Para solucionar el sistema planteado se eligió invertir la matriz 𝑀 y de allí
multiplicar 𝑀−1 con el vector 𝐴 para encontrar finalmente las variables del
vector .
5. Organizar y re-establecer las variables de flujo:
El desarrollo del proyecto en su mayoría fue implementado usando la
programación estructurada ya que esta permite, agilizar el desarrollo,
facilitar la lectura y comprensión de código y el cumplimiento de
estándares básicos de la codificación, por lo cual, es preferible seguir
este lineamiento mientras sea posible, sin embargo, es necesario
optimizar la implementación sobre el hardware, y para llevar esto a cabo
este paso es el encargado de solucionar los inconvenientes conceptuales
que se pudiesen presentar sobre la implementación.
Para hablar de optimización de la implementación en hardware es
necesario trabajar las variables y los accesos a ellas de la forma más
eficiente posible. Por ello este paso se encarga de traducir las variables
48
que se encuentran sobre posiciones de memoria estructuradas hacia
variables de acceso más veloz como simples vectores o matrices.
6. Iniciar configuración del Hardware
Para cumplir con el objetivo de filtrar señales mediante la reconstrucción
del sistema de filtrado es necesario configurar sobre la DSP todos los
elementos involucrados en esta tarea.
7. Filtrar mediante optimización del hardware
Para realizar una optimización del uso del hardware entre muchas
opciones resaltan las siguientes dos:
- Optimización mediantes hilos de procesamiento:
Consiste en implementar tareas y sincronizar resultados para que
estos procesos se lleven a cabo concurrentemente [10].
- Optimización de tiempo de ejecución realizada por el compilador:
Consiste en una optimización de código fuente realizada por el
compilador cuando se está traduciendo de C/C++ hacia
Assembler. Esto se realiza mediante el uso de los pragmas e
instrucciones para el compilador para que la traducción de código
que se realice resulte con la configuración deseada [20].
El método elegido para la codificación de este proyecto es la
optimización por compilador puesto que la arquitectura de la ADSP BF-
533 no permite el manejo de arreglos dinámicos de hilos.
Cabe hacer hincapié en que los algoritmos de filtrado se deben
implementar con los pragmas e instrucciones de compilador pertinentes
para que la traducción del código de C/C++ hacia Assembler logre la
optimización mencionada.
4.5 Implementación sobre Visual DSP++
En general la codificación de todo el proyecto sobre Visual DSP++ se realizó en
lenguaje C++ siguiendo los lineamientos de la programación estructurada por lo
cual el tratamiento de los diferentes procesos se pudieron encapsular en
estructuras y redefinición de operadores para facilitar la codificación, lectura del
código y el cumplimiento de estándares básicos de codificación.
49
El código fuente de todas las estructuras diseñadas y utilizadas en este
proyecto se muestra en el Anexo A.
A continuación se describirá cada uno de los pasos ya mencionados enfocados
hacia la codificación sobre el Visual DSP++:
1. Inicialización de Variables Fundamentales
La inicialización de las variables en el Visual DSP++ se realizó teniendo
en cuenta la necesidad de que ciertas variables fuesen estructuradas y
otras por optimización fueran simples vectores y matrices. A continuación
se mencionan algunas variables y sus respectivos tipos de inicialización:
- 4 Variables para el manejo de los 4 canales de entrada. Variables
inicializadas como float con valor por defecto 0.
- 4 Variables para el manejo de los 4 canales de salida. Variables
inicializadas como float con valor por defecto 0.
- Registro de configuración del códec 1836. Variable
sCodec1836TxRegs inicializada como short [] con los valores
para la configuración del códec de audio 1836.
- Buffer de tres posiciones para las muestras de entrada, variable
BE inicializada como arreglo de dos posiciones de datos tipo
flotante (float []).
- Matriz de Buferes de M x 2 posiciones para las muestras de
salida de cada subsistema. Variable MBS inicializada como
matriz de tipo flotante (*float []).
- Matriz de coeficientes de 2M x 3 para almacenar los coeficientes
de los subsistemas resultantes. Variable matrizCoef inicializada
como matriz de tipo flotante (*float []).
- Polinomio de numerador del sistema original a implementar.
Variable numSis implementada sobre la estructura Polinomio.
- Polinomio de denominador del sistema original a implementar.
Variable denSis implementada sobre la estructura Polinomio.
2. Lectura del sistema a Implementar:
50
Para la lectura del sistema a implementar se deben realizar el ingreso de
los coeficientes y su carga posterior a la memoria de la DSP como se
muestra abajo:
Para el caso de este trabajo estos coeficientes fueron ingresados
directamente sobre código.
3. Hallar los polos del sistema a implementar:
Para hallar los polos del sistema se puede seguir cualquier algoritmo que
pueda hallar las raíces reales y complejas de un polinomio con
coeficientes reales. En el caso de la implementación hecha sobre
VisualDSP++ se realizó mediante una variación para polinomios de
coeficientes reales del algoritmo general de Jenkins–Traub [21] publicado
por primera vez en 1970, este algoritmo es un método iterativo de
convergencia rápida y global. Cabe resaltar que su uso fue solo funcional
puesto que su estudio excede el propósito de este trabajo. El código
nativo de C++ utilizado se puede observar en el Anexo B, este tuvo que
ser adaptado para su correcto funcionamiento sobre la DSP sin embargo,
su flujo no fue alterado.
El diagrama de flujo mostrado a continuación describe el procedimiento
realizado para la implementación del algoritmo de Jenkins–Traub con el
fin de hallar los polos del sistema a implementar.
51
Es necesario llevar a cabo los casteos de las variables puesto que el
algoritmo de Jenkins-Traub es invocado como se encontró y no sigue el
patrón de codificación que se implementó en el proyecto (programación
estructurada).
4. Plantear y resolver el sistema de ecuaciones:
El proceso realizado para el planteamiento y solución del sistema de
ecuaciones para aplicar fracciones parciales al sistema original, es uno
de los núcleos fundamentales que posee este trabajo. Los numeradores
de los subsistemas se considerarán como polinomios de orden uno por lo
cual para hallar el valor de los términos de la columnas i e i+1 de la matriz
M, es necesario llevar a cabo la multiplicación de todos los polinomios
encontrados que no contengan las incógnitas de dichas columnas, esto
se realiza en la inicialización y el llenado de la variable mulPol la cual
almacena en su posición i-esima la multiplicación de todos los polinomios
diferentes al i-esimo polinomio. Todos los polinomios pertenecientes a los
polos hallados son almacenados en la variable polden.
A continuación se muestra una representación gráfica del algoritmo de
fracciones parciales implementado, es necesario tener en cuenta el
significado de cada una de las variables usadas en el esquema mostrado
por ello se describen a continuación:
52
- polden: Variable con todos los polinomios generados a partir del
denominador del sistema original.
- M: Matriz de relaciones.
- Mi: Matriz inversa a M.
- L: Vector de incógnitas.
- A: Vector con los coeficientes del numerador original (a excepción
del coeficiente del termino de mayor exponente, puesto que este
es la constante espectral).
- index: Índice de la fila de ceros en la matriz M.
- index2: Índice de la columna por la cual se cambió la columna 0.
- result: Variable resultado de tipo Fracciones.
53
Se debe especificar que la columna i e i+1 corresponden al coeficiente de
orden uno y de orden cero respectivamente pertenecientes al subsistema
j por lo cual los términos de las columnas i e i+1 son idénticos pero sus
posiciones difieren en una unidad en su índice de fila. Específicamente,
para la columna i+1 el valor de la fila j+1 es igual al valor de la fila j en la
columna i. Por lo anterior se ve que los términos de la matriz M son los
coeficientes de los polinomios contenidos en mulPol teniendo en cuenta
la regla de la diferencia en las filas para las columnas i e i+1. Abajo se
muestra la estructura general de la matriz M a partir de las
consideraciones mencionadas:
En algunas ocasiones al hallar los términos de la matriz M, su primera fila
contenía solo ceros y para corregir este inconveniente se eliminaba dicha
54
fila y una columna perteneciente a un subsistema de orden uno donde se
cumpliese que su termino 𝑀1𝑖! = 0. Se debe tener cuenta este hecho al
momento de hallar los valores de las incógnitas para corregir en ellas, la
modificación sobre la matriz M.
Para la inversión de la matriz M se empleó un método iterativo para el
cálculo de la inversa mediante los cofactores y el determinante de la
matriz, en varias ocasiones se presentó divergencia en dicho cálculo
cuando 𝑀00 era igual a cero. Para solucionar este inconveniente se debe
intercambiar la columna cero por la columna iesima donde 𝑀0𝑖! = 0, esto
produce que la incógnita de índice cero tenga el valor de la incógnita
iesima por lo cual se debe corregir esta situación al momento del cálculo
final de las incógnitas.
Finalmente para hallar el valor de las incógnitas se debe realizar la
multiplicación de la inversa por el vector de coeficientes del numerador y
corregir los posibles cambios realizados sobre la matriz M.
5. Organizar y re-establecer variables de flujo:
Las incógnitas halladas en el punto anterior son entregadas en una
estructura Fracciones por lo cual se deben cambiar las siguientes
variables a matrices y vectores:
- Matriz de coeficientes de los subsistemas.
- Arreglo con los órdenes de los subsistemas.
A continuación se describe la estructura de la matriz de coeficientes:
𝑀𝐶𝑂 =
[ 𝐶00
𝐷00
⋮
𝐶01
𝐷01
⋮
𝐶02
𝐷02
⋮𝐶𝑖0
𝐷𝑖0
⋮
𝐶𝑖1
𝐷𝑖1
⋮
𝐶𝑖2
𝐷𝑖2
⋮𝐶𝑀0
𝐷𝑀0
𝐶𝑀1
𝐷𝑀1
𝐶𝑀2
𝐷𝑀2]
(20)
𝐶𝑜𝑛 𝐶𝑖𝑗 𝑐𝑜𝑚𝑜 𝑒𝑙 𝑐𝑜𝑒𝑓𝑖𝑐𝑖𝑒𝑛𝑡𝑒 𝑑𝑒 𝑋2−𝑗 𝑑𝑒𝑙 𝑛𝑢𝑚𝑒𝑟𝑎𝑑𝑜𝑟 𝑑𝑒𝑙 𝑠𝑖𝑠𝑡𝑒𝑚𝑎 𝑖𝑒𝑠𝑖𝑚𝑜 𝑦
𝐷𝑖𝑗 𝑐𝑜𝑚𝑜 𝑒𝑙 𝑐𝑜𝑒𝑓𝑖𝑐𝑖𝑒𝑛𝑡𝑒 𝑑𝑒 𝑋2−𝑗 𝑑𝑒𝑙 𝑑𝑒𝑛𝑜𝑚𝑖𝑛𝑎𝑑𝑜𝑟 𝑑𝑒𝑙 𝑠𝑖𝑠𝑡𝑒𝑚𝑎 𝑖𝑒𝑠𝑖𝑚𝑜.
55
La variable que contiene los órdenes de los sistemas es un arreglo de
enteros donde la posición i-esima almacena el grado del subsistema i-
esimo.
6. Iniciar configuración del Hardware
La inicialización del hardware a utilizar sobre la DSP debe de tener en
cuenta los siguientes ítems a configurar:
- EBIU
- Códec 1836
- Puerto SPORT0
- DMA
- Configuración de Interrupciones
En el Anexo C se puede ver el código fuente de estas configuraciones
sobre la ADSP BF533 en el VisualDSP++.
7. Filtrar mediante optimización del hardware
Para filtrar optimizando el uso del hardware es necesario codificar
dándole las instrucciones adecuadas al compilador mediante el uso de
pragmas y de configuraciones de compilación.
A continuación se listan los pragmas utilizados:
- #pragma optimize_for_speed
- #pragma different_banks
- #pragma no_alias
Abajo las instrucciones dadas al compilador (se pasan como parámetros
al ejecutar el compilador):
- -O[0|1] seteada como -O1
- -extern-inline
- -Oa
La descripción de cada una de estas configuraciones se encuentra en la
referencia [20].
El código fuente para el filtrado de la señal puede verse a continuación
donde iChannel0LeftIn es la muestra actual, NUM_TAPS = 3,
iChannel0LeftOut variable de la muestra de salida y resultado variable
temporal para la sumatoria de las muestras de los subsistemas.
56
BE [NUM_TAPS-3]= BE [NUM_TAPS-2];
BE [NUM_TAPS-2]= BE [NUM_TAPS-1];
float resultado = 0;
#pragma optimize_for_speed
#pragma different_banks
#pragma no_alias
for( int sis= 0 ; sis < M ;sis++)
int numsis = 2*sis; //índice del numerador
int densis = numsis +1; //índice del denominador
float xx = MCO [numsis][0];
float yy = MCO [numsis][1];
float zz = MCO [numsis][2];
float vv = MCO [densis][1];
float ww = MCO [densis][2];
if(grado_sistemas[sis] == 2)
MBS [sis][2] =
(xx * iChannel0LeftIn) + (yy * BE [NUM_TAPS-2]) + (zz *
BE[NUM_TAPS-3]) - (vv * MBS [sis][1]) - (ww * MBS [sis][0]);
MBS [sis][0] = MBS [sis][1];
MBS [sis][1] = MBS [sis][2];
else
MBS [sis][2] =
(yy * iChannel0LeftIn + (zz * BE [NUM_TAPS-2]) - (ww *
MBS[sis][1]);
MBS [sis][1] = MBS [sis][2];
resultado = resultado + MBS [sis][2];
resultado = resultado + (BE [2] * CE);
iChannel0LeftOut = resultado;
4.6 Consecuencias Prácticas
En la Figura 21. Estructura genérica del banco de filtros para la reconstrucción del
sistema original se puede observar un banco de filtros donde su diseño es idéntico
al utilizado para la reconstrucción que se realiza de cualquier sistema. Obsérvese
que la entrada para todos los subsistemas es común en la fase de análisis. Luego
de llevarse a cabo el procesamiento de la muestra empieza la fase de síntesis del
banco de filtros y es importante resaltar que las entradas de esta fase son
57
independientes unas de otras.
Figura 21. Estructura genérica del banco de filtros para la reconstrucción del sistema original.
Cuando se desea implementar un sistema de filtrado digital de orden N normalmente se deben tener dos buffers de N posiciones para almacenar las muestras de salida y de entrada que se van obteniendo al procesar cada muestra.
La metodología planteada presenta un mayor rendimiento en memoria puesto
que al tener solo sistemas de orden uno y dos se tiene que los buferes de entrada y salida serán solo de una o dos posiciones respectivamente, además, como las entradas de todos los subsistemas deben ser las mismas entonces solo es necesario tener un único buffer de entrada de dos posiciones, y buffers de salida de una o dos posiciones únicos para cada subsistema. La Tabla 3 muestra el análisis de uso de memoria de la metodología tradicional y de la que se plantea en este trabajo.
Tabla 3. Cálculo del uso de memoria de la implementación por medio de la metodología tradicional y la planteada.
Uso de memoria para un sistema de orden N.
Metodología Utilizada
Posiciones del Buffer Entrada
Posiciones del Buffer Salida
Total Posiciones de Memoria
Tradicional N N 2N
Nueva 2 N N+2
58
5. RESULTADOS
Los resultados expuestos aquí resaltan las ventajas y los problemas
encontrados con la implementación de la metodología sobre el procesador
digital de señales ADSP BF-533.
Consideraciones a tener en cuenta al examinar los resultados presentados:
- Las respuestas en frecuencia teóricas fueron halladas mediante el
software Matlab®.
- Algunos de los filtros analizados son comúnmente utilizados en la
literatura puesto que se desea evaluar la pertinencia de la
metodología desarrollada.
- Las respuestas en frecuencia halladas en la práctica fueron
medidas mediante la tarjeta de adquisición NI MyDAQ utilizando
la función pre-programada BODE descrita anteriormente.
- Puesto que los valores de las muestras, que entran y salen del
códec 1836 de la DSP son re-escalados por el hardware para
evitar problemas de saturación, los valores máximos de ganancia
en la banda pasante de la respuesta en frecuencia para la
magnitud hallada en la práctica y teóricamente no coincidirán, sin
embargo, la forma de la respuesta en frecuencia si deberán ser
idénticas.
- Se tomarán dos frecuencias y sus respectivas magnitudes para el
análisis del delta de magnitud en las respuestas en frecuencia
teórica y práctica. Para escoger estas frecuencias se deberá tener
en cuenta solo las frecuencias que no se vean afectadas por
parámetros de hardware.
- La mínima ganancia que la NI MyDAQ lee es de
aproximadamente -55 dB.
59
Resultado 1.
La primera prueba se realizó contra un filtro Butterworth que se ha
implementado múltiples veces en la Universidad Tecnológica de Pereira cuando
se inicia el estudio de la ADSP BF-533, se trata del ejemplo mostrado en el anexo
F de [22]. La función de transferencia discreta de este filtro se muestra a
continuación:
𝐻(𝑧) =0.0002087 ∙ 𝑍4 + 0.0008349 ∙ 𝑍3 + 0.001252 ∙ 𝑍2 + 0.0008349 ∙ 𝑍1 + 0.0002087
𝑍4 + 3.321 ∙ 𝑍3 − 4.183 ∙ 𝑍2 + 2.365 ∙ 𝑍 − 0.5053(21)
Descripción
del filtro: Filtro pasa bajos de 4to orden Butterworth con frecuencia de corte de 2KHz
Velocidad de
Caída: 80 dB/década Frecuencia de muestreo: 98KHz
Respuesta
en
Frecuencia
para la
magnitud
hallada
teóricamente
Respuesta
en
Frecuencia
para la
magnitud
medida en la
práctica
Frecuencia 1: 4000 Hz Frecuencia 2: 8000 Hz
Delta de ganancia teórico: 24.7 dB Delta de ganancia hallado
en la práctica: 21.4 dB
Ganancia en banda pasante
hallada en la práctica: ≈ -10.98 dB
Frecuencia de corte
hallada en la práctica: ≈ 1995.26 Hz
Error Porcentual sobre el
delta de ganancia: ≈ 13%
Error porcentual sobre la
frecuencia de corte: ≈ 0.23%
60
Resultado 2.
El desarrollo de aplicaciones para la industria es uno de los lineamientos
que toda ingeniería debe considerar. Aquí se presenta un filtro diseñado e
implementado por Behringer® usado en diversos modelos de consolas
mezcladores de audio profesional que esta marca ofrece [23]. En la ecuación
22 se puede observar la función de transferencia del filtro implementado.
𝐻(𝑧) =0.49981 ∙ 𝑍4 − 1.999 ∙ 𝑍3 + 2.9988 ∙ 𝑍2−1.9992 ∙ 𝑍1 + 0.4998149
𝑍4 − 2.63862 ∙ 𝑍3 + 2.76930 ∙ 𝑍2 − 1.3392 ∙ 𝑍 + 0.24982(22)
Descripción
del filtro: Filtro pasa altos de 4to orden Butterworth con frecuencia de corte de 8KHz
Velocidad de
Caída: 80 dB/década Frecuencia de muestreo: 98KHz
Respuesta
en
Frecuencia
para la
magnitud
hallada
teóricamente
Respuesta
en
Frecuencia
para la
magnitud
medida en la
práctica
Frecuencia 1: 4000 Hz Frecuencia 2: 8000 Hz
Delta de ganancia teórico: 21.83 dB Delta de ganancia hallado
en la práctica: 15.23 dB
Ganancia en banda pasante
hallada en la práctica: ≈ -9.87 dB
Frecuencia de corte
hallada en la práctica: ≈ 7244.3 Hz
Error Porcentual sobre el
delta de ganancia: ≈ 30.23%
Error porcentual sobre la
frecuencia de corte: ≈ 9.44%
61
Resultado 3.
En la actualidad el estudio de la electrofisiología del cuerpo humano es un área
de gran interés y una de sus mayores herramientas son las señales ECG. Lin Xu,
Chiara Rabotti y Massimo Mischi miembros de la IEEE en 2012 usaron en [24]
diferentes tipos de filtros para acoplar y poder procesar correctamente señales
ECG. A continuación se muestra uno de los filtros utilizados en para discriminar
las frecuencias de interés sobre el ruido en las señales de interés. La función de
transferencia implementada es mostrada a continuación:
𝐻(𝑧) =
1.463 ∙ 10−7 ∙ 𝑍5 + 7.316 ∙ 10−7 ∙ 𝑍4 + 1.463 ∙ 10−6 ∙ 𝑍3 + 1.463 ∙ 10−6 ∙ 𝑍2+7.316 ∙ 10−7 ∙ 𝑍1 + 1.463 ∙ 10−7
𝑍5 − 4.714 ∙ 𝑍4 + 8.8968 ∙ 𝑍3 − 8.4024 ∙ 𝑍2 + 3.9709 ∙ 𝑍 − 0.7512(23)
Descripción
del filtro: Filtro pasa bajos de 5to orden Butterworth con frecuencia de corte de 450 Hz
Velocidad de
Caída: 100 dB/década Frecuencia de muestreo: 32KHz
Respuesta
en
Frecuencia
para la
magnitud
hallada
teóricamente
Respuesta
en
Frecuencia
para la
magnitud
medida en la
práctica
Frecuencia 1: 550 Hz Frecuencia 2: 910 Hz
Delta de ganancia teórico: 21.42 dB Delta de ganancia hallado
en la práctica: 23.32 dB
Ganancia en banda pasante
hallada en la práctica: ≈ -8.15 dB
Frecuencia de corte
hallada en la práctica: ≈ 422.8 Hz
Error Porcentual sobre el
delta de ganancia: ≈ 8.87%
Error porcentual sobre la
frecuencia de corte: ≈ 6.04%
62
Resultado 4.
Las aplicaciones a la ingeniería de sonido resultantes de la unión de
estudios interdisciplinares como el desarrollo de software y la física de las
ondas mecánicas proveen grandes aportes en la actualidad para la innovación
de los dispositivos tecnológicos usados cotidianamente. Jinquan Zhang y
Hongxia Wang de la Escuela de Ciencias de la información y Tecnología
utilizaron diferentes filtros para estudiar sus efectos sobre el proceso de
WaterMarking2 en archivos de audio [25]. A continuación se muestra uno de
estos filtros con su respectiva función de transferencia implementado mediante
la metodología propuesta.
𝐻(𝑧) =
1.105 ∙ 10−5 ∙ 𝑍10 + 0.0001 ∙ 𝑍9 + 0.0004 ∙ 𝑍8 + 0.0013 ∙ 𝑍7+0.0023 ∙ 𝑍6+0.0027 ∙ 𝑍5 …
𝑍10 − 4.986 ∙ 𝑍9 + 11.936 ∙ 𝑍8 − 17.742 ∙ 𝑍7 + 17.973 ∙ 𝑍6 − 12.886 ∙ 𝑍5 …
…0.0023 ∙ 𝑍4 + 0.0013 ∙ 𝑍3 + 0.0004 ∙ 𝑍2 + 0.0001 ∙ 𝑍 + 1.105 ∙ 10−5
…+ 6.5932 ∙ 𝑍4 − 2.3690 ∙ 𝑍3 + 0.5706 ∙ 𝑍2 − 0.0830 ∙ 𝑍 + 0.0055(24)
Descripción
del filtro:
Filtro pasa bajos de 10mo orden
Butterworth con frecuencia de corte de 6KHz
Velocidad de
caida: -200 dB Frecuencia de muestreo: 48KHz
Respuesta en
Frecuencia
para la
magnitud
hallada
teóricamente
2 Proceso de estampado digital con el fin de poner de manifiesto el uso ilegal de archivos digitales.
63
Respuesta en
Frecuencia
para la
magnitud
medida en la
práctica
Frecuencia 1: 6 KHz Frecuencia 2: 8 KHz
Delta de ganancia teórico: 50.89 dB Delta de ganancia
hallado en la práctica: 30.93 dB
Ganancia en banda pasante
hallada en la práctica: ≈ -8.71 dB
Frecuencia de corte
hallada en la práctica: ≈ 5703.12 Hz
Error Porcentual sobre el
delta de ganancia: ≈ 39.22 %
Error porcentual sobre la
frecuencia de corte: ≈ 4.94 %
Resultado 5.
Abajo se detallan cuatro diseños propios de filtros para poner a prueba la
metodología y con el fin de mostrar al lector el comportamiento de los diferentes
tipos de filtros implementados mediante la metodología.
𝐻(𝑧) =
0.4956 ∙ 𝑍6 − 2.855 ∙ 𝑍5 + 6.970 ∙ 𝑍4 − 9.220 ∙ 𝑍3 + 6.9700 ∙ 𝑍2−2.8553 ∙ 𝑍1 + 0.4956
𝑍6 − 4.4613 ∙ 𝑍5 + 8.4376 ∙ 𝑍4 − 8.7329 ∙ 𝑍3 + 5.2486 ∙ 𝑍2 − 1.7368 ∙ 𝑍 − 0.2451 (25)
Descripción
del filtro:
Filtro rechaza banda de 6to orden
Butterworth con frecuencias de corte de 500 Hz y 4KHz
Velocidad de
Caída: 60 dB/década Frecuencia de muestreo: 32KHz
Respuesta
en
Frecuencia
para la
magnitud
hallada
teóricamente
64
Respuesta
en
Frecuencia
para la
magnitud
medida en la
práctica
Frecuencia 1: 500 Hz Frecuencia 2: 860 Hz
Delta de ganancia teórico: 19.19 dB Delta de ganancia hallado
en la práctica: 21.25 dB
Ganancia en banda pasante
hallada en la práctica: ≈ -8.08 dB
Frecuencia de corte
hallada en la práctica: ≈ 501.19 Hz
Error Porcentual sobre el
delta de ganancia: ≈ 10.7%
Error porcentual sobre la
frecuencia de corte: ≈ 0.23%
𝐻(𝑧) =
0.014 ∙ 𝑍8 − 0.071 ∙ 𝑍7 + 0.178 ∙ 𝑍6 − 0.29 ∙ 𝑍5 + 0.338 ∙ 𝑍4 − 0.29 ∙ 𝑍3 + 0.178 ∙ 𝑍2−0.071 ∙ 𝑍1 + 0.014
𝑍8 − 6.465 ∙ 𝑍7 + 19.1 ∙ 𝑍6 − 33.51 ∙ 𝑍5 + 38.03 ∙ 𝑍4 − 28.56 ∙ 𝑍3 + 13.85 ∙ 𝑍2 − 3.964 ∙ 𝑍 + 0.5135 (26)
Descripción
del filtro:
Filtro pasa bajos de 8vo orden
elíptico con frecuencia de corte de 3KHz
Ganancia en
banda rechazo: -40 dB Frecuencia de muestreo: 32KHz
Respuesta en
Frecuencia
para la
magnitud
hallada
teóricamente
65
Respuesta en
Frecuencia
para la
magnitud
medida en la
práctica
Frecuencia 1: 3 KHz Frecuencia 2: 3.14 KHz
Delta de ganancia teórico: 39.6 dB Delta de ganancia
hallado en la práctica: 37.28 dB
Ganancia en banda pasante
hallada en la práctica: ≈ -10.34 dB
Frecuencia de corte
hallada en la práctica: ≈ 2985.4 Hz
Error Porcentual sobre el
delta de ganancia: ≈ 5.85%
Error porcentual sobre la
frecuencia de corte: ≈ 0.48%
𝐻(𝑧) =
0.000194 ∙ 𝑍8 − 0.00077 ∙ 𝑍6 + 0.00116 ∙ 𝑍4 − 0.000776 ∙ 𝑍2 + 0.000194
𝑍8 − 7.1357 ∙ 𝑍7 + 22.47 ∙ 𝑍6 − 40.79 ∙ 𝑍5 + 46.69 ∙ 𝑍4 − 34.50 ∙ 𝑍3 + 16.07 ∙ 𝑍2 − 4.318 ∙ 𝑍 + 0.5121 (27)
Descripción
del filtro:
Filtro pasa banda de 8to orden
Butterworth con frecuencias de corte de 700 Hz y 2KHz
Velocidad de
Caída: 80 dB/década Frecuencia de muestreo: 32KHz
Respuesta
en
Frecuencia
para la
magnitud
hallada
teóricamente
66
Respuesta
en
Frecuencia
para la
magnitud
medida en la
práctica
Frecuencia 1 (Low): 2.5 KHz Frecuencia 1 (High): 430 Hz
Frecuencia 2 (Low): 4 KHz Frecuencia 2 (High): 700 Hz
Delta de ganancia teórico
(Low): 23 dB
Delta de ganancia teórico
(High): 23.89 dB
Delta de ganancia hallado en
la práctica (Low): 23.11 dB
Delta de ganancia hallado
en la práctica (High): 22.5 dB
Frecuencia de corte hallada
en la práctica (Low): ≈ 1925.6 Hz
Frecuencia de corte hallada
en la práctica (High): ≈ 637.3 Hz
Error Porcentual sobre el
delta de ganancia (Low): ≈ 0.47%
Error porcentual sobre el
delta de ganancia (High): ≈ 5.81%
Error porcentual sobre la
frecuencia de corte (Low): ≈ 3.72%
Error porcentual sobre la
frecuencia de corte (High): ≈ 8.95%
𝐻(𝑧) =
10−5 ∙ 𝑍8 + 0.0001 ∙ 𝑍7 + 0.0004 ∙ 𝑍6 + 0.0008 ∙ 𝑍5 + 10−3 ∙ 𝑍4 + 0.0008 ∙ 𝑍3 + 0.0004 ∙ 𝑍2+0.0001 ∙ 𝑍1 + 10−5
𝑍8 − 4.984 ∙ 𝑍7 + 11.261 ∙ 𝑍6 − 14.94 ∙ 𝑍5 + 12.68 ∙ 𝑍4 − 7.034 ∙ 𝑍3 + 2.480 ∙ 𝑍2 − 0.507 ∙ 𝑍 + 0.046 (28)
Descripción
del filtro:
Filtro pasa bajos de 8vo orden
Butterworth con frecuencia de corte de 3KHz
Velocidad de
Caída: 160 dB/década Frecuencia de muestreo: 32KHz
Respuesta
en
Frecuencia
para la
magnitud
hallada
teóricamente
67
Respuesta
en
Frecuencia
para la
magnitud
medida en la
práctica
Frecuencia 1: 3.5 KHz Frecuencia 2: 5 KHz
Delta de ganancia teórico: 27.6 dB Delta de ganancia hallado
en la práctica: 25.33 dB
Ganancia en banda pasante
hallada en la práctica: ≈ -10.32 dB
Frecuencia de corte
hallada en la práctica: ≈ 2642.69 Hz
Error Porcentual sobre el
delta de ganancia: ≈ 8.22%
Error porcentual sobre la
frecuencia de corte: ≈ 11.9%
Resultado 6.
El diseño de filtros digitales debe tener en cuenta factores prácticos como lo
son: la frecuencia de muestreo, reloj principal del procesador, restricciones de
hardware, diversos tipos de ruidos, etc. El siguiente resultado muestra un mal
diseño adrede, con el fin de ilustrar las restricciones de hardware que la ADSP BF-
533 presenta. La ecuación 29 muestra la función de transferencia del filtro
implementado.
𝐻(𝑧) =
4.60 ∙ 10−8 ∙ 𝑍8 + 3.68 ∙ 10−7 ∙ 𝑍7+1.28 ∙ 10−6 ∙ 𝑍6+2.57 ∙ 10−6 ∙ 𝑍5 …
𝑍8 − 6.65 ∙ 𝑍7 + 19.49 ∙ 𝑍6 − 32.759 ∙ 𝑍5 …
…3.22 ∙ 10−6 ∙ 𝑍4 + 2.57 ∙ 10−6 ∙ 𝑍3 + 1.28 ∙ 10−6 ∙ 𝑍2 + 3.68 ∙ 10−7 ∙ 𝑍 + 4.60 ∙ 10−8
…+ 34.55 ∙ 𝑍4 − 23.42 ∙ 𝑍3 + 9.961 ∙ 𝑍2 − 2.429 ∙ 𝑍 + 0.26003 (29)
68
Descripción
del filtro:
Filtro pasa bajos de 8vo orden
Butterworth con frecuencia de corte de 4KHz
Ganancia en
banda rechazo: -40 dB Frecuencia de muestreo: 96KHz
Respuesta en
Frecuencia
para la
magnitud
hallada
teóricamente
Respuesta en
Frecuencia
para la
magnitud
medida en la
práctica
Frecuencia de corte hallada
en la práctica: ≈ 2293.55 Hz
Error porcentual sobre
frecuencia de corte ≈ 42.66 %
69
La causa del mal funcionamiento del filtro implementado tiene varios factores a
considerar:
- Frecuencia de muestreo:
Puesto que los métodos de diseño tradicionales de filtros digitales emplean el
rango de posibles frecuencias entre 0 y 1 (donde 1 es la mitad de la frecuencia de
muestreo) los coeficientes del numerador de la función de transferencia de un filtro
-cuya frecuencia de corte no se desee alterar- cambiarán en función de la
frecuencia de muestreo, y en general, se harán más pequeños cada vez que se
aumente la frecuencia de muestreo conservando la frecuencia del filtro inalterada.
- Orden elevado del filtro:
A medida que la selectividad frecuencial de un filtro digital aumenta, los
coeficientes en el numerador de este tenderán a cero por lo cual, a órdenes muy
elevados, es necesario contar con el hardware adecuado para el procesamiento
del filtro.
- Cifras significativas de los float’s en la ADSP BF-533:
En la implementación sobre la ADSP BF-533 las variables con las cuales se
hallaban los coeficientes de los subsistemas eran de tipo float y por esta razón
cuentan con solo ocho cifras significativas a la hora de su uso.
El sistema de este resultado en particular, en el numerador de su función de
transferencia, contaba con coeficientes del orden de 10−8 por lo cual, al hallar su
equivalente mediante los subsistemas se encontraban aproximaciones erróneas lo
cual genera el mal resultado del filtro. Estas malas aproximaciones son debidas a
la precisión con que las variables float se manejan sobre el hardware.
70
6. CONCLUSIONES
Las aplicaciones de ingeniera actualmente tienen embebidas consigo
diversos tipos de filtros para rechazar componentes frecuenciales indeseados
sobre las señales de interés. La metodología presentada demostró mediante
los resultados obtenidos que ofrece una forma de implementar filtros LTI de
manera sencilla, rápida, eficaz, en tiempo real y sobre un sistema embebido
mediante el uso óptimo del hardware disponible para el procesamiento de
cualquier aplicación ingenieril. Gracias a la eficiente realización que la
metodología lleva a cabo, las aplicaciones que hagan uso de ella podrán tener
más capacidad de procesamiento para las demás tareas que implementando
los filtros normalmente, puesto que automáticamente, y sin ninguna
intervención del usuario, se optimiza la realización de filtro deseado.
El uso de la metodología presenta grandes ventajas para la construcción de
la respuesta en frecuencia deseada ya que, al reconstruir el sistema
inicialmente diseñado a partir de subsistemas, permite el manejo individual de
cada uno de ellos, dándole al usuario la capacidad de seleccionar ganancias
particulares a cada subsistema y así poder construir una respuesta en
frecuencia optimizada y particularmente parametrizada en función de la
aplicación.
La implementación de la metodología presenta problemas al momento de
reconstruir sistemas que cuenten con polos anidados, esto es producto del uso
del algoritmo de fracciones parciales en la metodología, puesto que la
reconstrucción de sistemas con esta característica podría conllevar la
implementación de filtros de orden mayor a dos, lo cual no está concebido en la
codificación realizada. Formalizando se tiene que los sistemas que cumplan
con la siguiente expresión no podrán reconstruirse mediante el análisis
presentado:
H(z) =A(z)
(z ± a)n ∗ ((z + bi)(z − bi))m ∗ B(z) con n > 2,m > 2
En el transcurso de este trabajo se presentaron diversos problemas
referentes al hardware utilizado, el mayor ellos, fue que las posibles
frecuencias de muestreo que la ADSP BF-533 posee son: 32KHz, 44KHz,
48KHz y 96KHz. A raíz de lo anterior, la implementación de filtros que fuesen
muy selectivos frecuencialmente solo se podían reconstruir si los coeficientes
71
de su función de transferencia en el numerador eran de órdenes mayores a 10−8
producto de la precisión de las variables utilizadas. Por este motivo al momento de
hacer uso de la metodología sobre cualquier dispositivo con el fin de filtrar señales
se debe tener en cuenta los siguientes factores:
Frecuencia de Muestreo del ADC:
Experimentalmente se encontró que la frecuencia de muestreo con la cual
la metodología funciona bien debe ser como máximo 70 veces mayor a la
frecuencia de corte del filtro a implementar.
Cifras significativas que las variables usadas posean:
Al momento de hallar los subsistemas que reconstruirán el filtro, la precisión
de las variables allí utilizadas es la que determina el valor mínimo que los
coeficientes podrán tener en los sistemas a implementar.
Orden de los coeficientes del numerador del sistema
El valor mínimo de los coeficientes en los numeradores del sistema a
implementar viene dado por el valor del delta entre los posibles valores del
tipo de variable utilizada en la codificación del proyecto.
Los problemas de implementación que presentó el hardware de la ADSP BF-
533 fueron incisos claves para determinar las limitantes de la implementación
metodológica y gracias a ellos se pudo delimitar cuales son los factores y
requisitos que más impactan la implementación de la metodología. Se debe
resaltar la necesidad de mejorar la metodología para la determinación y
autoconfiguración de la frecuencia de muestreo del sistema a implementar, ya que
esto podría mejorar de manera drástica el rango de sistemas que la metodología
puede implementar evadiendo algunas de las restricciones que el hardware
imponga.
A pesar de las restricciones impuestas por el hardware utilizado y los
problemas de codificación, se mostró mediante los resultados obtenidos que la
metodología cumple su propósito de manera eficaz en la reconstrucción de
sistemas de filtrado mediante la superposición de las respuestas en frecuencia de
los subsistemas hallados mediante la aplicación de fracciones parciales al sistema
de filtrado original.
72
7. BIBLIOGRAFÍA
[1] DAVID P. MORGAN, KEVIN W. LEARY. “A Dynamic Programming-Based Dur
Recognizer Implemented on Analog Device's DSP”, 1986. Brown University, One
Technology Way, Norwood.Analog Device, USA.
[2] KATSUHIKO HAYASHI, KAUSHAL K. DHAR, KAZUNORI SUGAHARA and
KOTARO HIRANO, “Design of High-Speed Digital Filters Suitable for Multi-DSP
Implementation”, Departaments of Engienering Electronic and Electric from Kobe
University, Japon. IEEE Transactions On Circuits And Systems, Vol 33, No 2,
Febrero 1986.
[3] WALTER DEMMER, JOACHIM SPEIDEL. “Programmable Digital Filter For
High Speed Aplications”. Philips Komunikations Industrie AG. CH2458 in
ISCAS’88.
[4] JONNE POIKONEN Y ARI PAASIO, “A Ranked Order Filter Implementation for
Parallel Analog Processing”, Departament of Information Technology of Tourku
University, Finland. IEEE Transactions On Circuits And Systems – Regular Papers,
Vol. 51, No 5, Mayo 2004.
[5] GUO GAIZHI, ZHANG PENGJU, YU ZONGZUO Y WANG HAILONG, “Design
and Implementation of FIR Digital Wave Filter Based on DSP”, College of
Computer and Information Technology, Inner Mongolia Normal University, Huhhot
China.978-14244-5824, IEEE 2010.
[6] ZHANG JIAN, PAN LING Y DING QUAN, “IIR Digital Filter Design Based on
DSP and SOPC”, School of Electric Power and Automation Engineering, Shanghai
University of Electric Power, China. International Symposium on Intelligence
Information Processing and Trusted Computing, DOI 10.1109, IEEE 2010.
[7] Code Composer Studio User’s Guide, Texas Instruments, Febrary 2000. Texas
Instruments.
[8]La Máquina de Turing, “http://www.sinewton.org/numeros/numeros/43-
44/Articulo33.pdf”, Manuel Fonseca de “Las matemáticas del siglo XX”.
[9]YAO LI Y XUEHUA ZHANG, “Parallel Evolutionary Method for Analog State
Variable Transconductance”, Physics School Beihua University From Jilin City,
China. International Conference On Mechatronic Science, Electric Engineering And
Computer, Jilin, China, August 2011.
73
[10] P. J. Dietel, Como Programar En Java, Pearson Prentice Hall, 7ma Ed, pp.
925-983.
[11] ALAN V. OPPENHEIM, ALAN S. WILLSKY, S. HAMID NAWAB, Señales y sistemas, Person Prentice Hall, 2da Ed, pp 651. [12] H. Nyquist, Certain topics in telegraph transmission theory, Trans. AIEE, vol. 47, pp 617-644, april 1928. [13] ALAN V. OPPENHEIM, ALAN S. WILLSKY, S. HAMID NAWAB, Señales y sistemas, Person Prentice Hall, 2da Ed, pp 33-35. [14] MURRAY R. SPIEGEL, Transformadas de Laplace, Serie Schaum, Mc Graw Hill, 1ra Ed, pp 1-6, 1997.
[15] L. MILIC, U. KAMSEK AND A. JOCIC, “Analysis off filter bank transfer functions”, Pupin Institute, Mihajlo Belgrade, IEE Proceedings, Vol. 134, pp 23-30, February 1987.
[16] ANALOG DEVICES INC. ®, Blackfin Embedded Processor ADSP-BF531/ ADSP-BF532/ ADSP-BF533, 1ra Ed, 2011.
[17] ANALOG DEVICES INC. ®, Blackfin Embedded Processor ADSP-BF531/ ADSP-BF532/ ADSP-BF533, 1ra Ed, pp 5, 2011.
[18] D. M. Carvajal Ramirez y P. F. Páramo Vaca, Desarrollo de una guía práctica de laboratorio par la materia de procesamiento digital de señales utilizando la tarjeta "ADZS-TS201S" de Analog Devices, Sangolquí - Ecuador: Escuela Politécnica del Ejercito, 2009.
[19] NATIONAL INSTRUMENTS ®, USER GUIDE NI MyDAQ, 2da Ed, 2014.
[20] VISUAL DSP++ 5.0, C/C++ Compiler and Library Manual for Blackfin ® Processors, Rev 5.4, enero 2011.
[21] M. A. JENKINS, J. F. TRAUB, “A three-stage variable-shift iteration for polynomial zeros and its relation to generalized Rayleigh iteration”, Washington University, vol. 14, pp 252-263, Febrero 1970.
[22] JUAN D. VALENCIA, “IMPLEMENTACIÓN DEL ALGORITMO DE PAN-TOMPKINS PARA LA DETECCIÓN DE COMPLEJOS QRS EN EL DSP ANALOG DEVICES ADSP-BF533”, Universidad Tecnológica de Pereira, pp 133-152, 2012.
[23] Behringer ®, “User Manual XENYX X1204USB / 1204USB”, 1ra Ed,North Creek Parway USA.
74
[24] Lin Xu, Chiara Rabotti, Massimo Mischi, “Novel Vibration-Exercise Instrument with Dedicated Adaptative Filtering for Electromyographic Investigation of Neuromuscular Activation”, Miembros de la IEEE, vol. 21, No 2, pp 275-282, Marso 2013.
[25] Jinquan Zhang, Hongxia Wang, “Audio Watermarking Scheme Resistant to Both Random Cropping and Low-Pass Filtering”, Escuela de Ciencias de la Información y Tencnología , vol. 21, No 2, pp 275-282, Marso 2013.
75
8. ANEXOS
Anexo A. Código fuente de las estructuras diseñadas e implementadas.
//ESTRUCTURAS //Estructura para numeros complejos struct Complex double Real; double Comp; bool Cor; Complex () Real=0; Comp=0; Cor=true; Complex (double r, double c) Real = r; Comp = c; Cor=true; ; //Estructura para el manejo de filas de matrices struct Fila int Count; bool Cor; double * Colm; Fila() Count=0; Colm = new double[1]; Cor=true; Colm[0]=0; Fila(int n)
Count=n; Colm = new double[n]; Cor=true; for(int i=0;i<n;i++) Colm[i]=0;
; //Estructura matriz donde esta contenida un arreglo de filas para cumplir con la //estructura de matriz struct Matriz int Fil; int Col; bool Cor; Fila *Filas; Matriz(int fil,int col) Fil=fil; Col=col; Filas = new Fila[fil]; for (int i =0; i<fil;i++) Filas[i] = Fila(col); Cor = true; ;
76
//Estructura para el manejo de raices de los polimonios struct Raices int Count; Complex * Raiz; bool Cor; Raices(int n) Count = n; Raiz = new Complex[n]; for(int i = 0; i<n;i++) Raiz[i]= Complex(); Cor = true; ; //Estructura para el manejo de polinomios. La posicion i es el coeficiente de //x^(Grado-i) donde Grado es el Grado del prolinomio de la instancia struct Polinomio int Grado; double * Coef; Polinomio() Grado = 0; Coef = new double[1]; Coef[0] = 1; Polinomio(int g) Grado = g; Coef = new double[g+1]; for(int i = 0;i<g+1;i++) Coef[i]=0; ; //Estructura para el manejo de funciones racionales struct Racional Polinomio Numer; Polinomio Denom; bool Cor; Racional () Numer = Polinomio(0); Denom = Polinomio(0); Denom.Coef[0]=1; Cor = true;
77
Racional (int num_n, int num_d) Numer = Polinomio(num_n); Denom = Polinomio(num_d); Numer.Coef[0]=0; Cor = true; ; //Estructura para el manejo de varias funciones racionales. Se utiliza para el //manejo del resultado de las fracciones parciales struct Fracciones Racional * Frac; double k; int Num; bool Cor; Fracciones() Num = 1; Cor = true; Frac = new Racional[1]; Frac[0] = Racional(); k=0; Fracciones(int num) Num = num; Cor = true; Frac = new Racional[num]; for(int i = 0;i< num ; i++) Frac[i] = Racional(); k=0; ; //Estructura auxiliar para el manejo de matrices, esta es la que se usa para hallar //la inversa de la matriz struct Matriz2 int actualsize; double *data; Matriz2(int n) actualsize=n; data = new double[n*n]; for (int i =0; i<n*n;i++)
78
data[i] = 0; ; //OPERADORES //Operador suma para estructuras de tipo Complex Complex operator + (Complex a, Complex b) return Complex(a.Real+b.Real, a.Comp+b.Comp); //Operador resta para estructuras de tipo Complex Complex operator - (Complex a, Complex b) return Complex(a.Real-b.Real, a.Comp-b.Comp); //Operador multiplicacion para estructuras de tipo Complex Complex operator * (Complex a, Complex b) return Complex( (a.Real*b.Real) - (a.Comp*b.Comp), (a.Real*b.Comp) + (a.Comp*b.Real) ); //Operador division para estructuras de tipo Complex Complex operator / (Complex a, Complex b) Complex Temp = Complex(); double doubletemp = (b.Real*b.Real) + (b.Comp*b.Comp); Temp.Real = ((a.Real*b.Real) + (a.Comp*b.Comp) )/doubletemp; Temp.Comp = ((a.Comp*b.Real) - (a.Real*b.Comp))/doubletemp; return Temp; //POLINOMIO //Operador multiplicacion para estructuras de tipo Polinomio Polinomio operator *(Polinomio a, Polinomio b) int i=0,j=0; Polinomio temp = Polinomio(a.Grado + b.Grado); if(a.Grado == 0) for(i=0;i<b.Grado+1;i++) temp.Coef[i]=b.Coef[i]*a.Coef[0]; else if(b.Grado == 0) for(i=0;i<a.Grado+1;i++) temp.Coef[i]=a.Coef[i]*b.Coef[0]; else for(i=0;i<a.Grado+1;i++) for(j=0;j<b.Grado+1;j++)
79
temp.Coef[i+j]=temp.Coef[i+j]+(a.Coef[i]*b.Coef[j]); int gr = Grado_Real_Pol(temp); Polinomio ress = Polinomio(gr); ress = Rectificar_Pol(temp,gr); return temp; //Operador suma para estructuras de tipo Polinomio Polinomio operator +(Polinomio a, Polinomio b) int orden,dif,i=0; Polinomio temp = Polinomio(0); if(a.Grado>b.Grado) orden = a.Grado; dif = orden-b.Grado; temp = Polinomio(orden); for(i=0;i<dif;i++) temp.Coef[i]=a.Coef[i]; for(i=dif;i<orden+1;i++) temp.Coef[i]=a.Coef[i]+b.Coef[i-dif]; else orden = b.Grado; dif = orden-a.Grado; temp = Polinomio(orden); for(i=0;i<dif;i++) temp.Coef[i]=b.Coef[i]; for(i=dif;i<orden+1;i++) temp.Coef[i]=a.Coef[i-dif]+b.Coef[i]; int gr = Grado_Real_Pol(temp); Polinomio ress = Polinomio(gr); ress = Rectificar_Pol(temp,gr); return ress; //Operador resta para estructuras de tipo Polinomio Polinomio operator -(Polinomio a, Polinomio b) int orden,dif,i=0; Polinomio temp = Polinomio(0);
80
if(a.Grado>b.Grado) orden = a.Grado; dif = orden-b.Grado; temp = Polinomio(orden); for(i=0;i<dif;i++) temp.Coef[i]=a.Coef[i]; for(i=dif;i<orden+1;i++) temp.Coef[i]=a.Coef[i]-b.Coef[i-dif]; else orden = b.Grado; dif = orden-a.Grado; temp = Polinomio(orden); for(i=0;i<dif;i++) temp.Coef[i]=-b.Coef[i]; for(i=dif;i<orden+1;i++) temp.Coef[i]=a.Coef[i-dif]-b.Coef[i]; int gr = Grado_Real_Pol(temp); Polinomio ress = Polinomio(gr); ress = Rectificar_Pol(temp,gr); return ress; //RACIONAL //Operador suma para estructuras de tipo Racional Racional operator +(Racional a,Racional b) int ordsup = 0; int ordinf = 0; if((a.Numer.Grado + b.Denom.Grado) >= (a.Denom.Grado + b.Numer.Grado)) ordsup = a.Numer.Grado + b.Denom.Grado; else ordsup = a.Denom.Grado + b.Numer.Grado; ordinf = a.Denom.Grado + b.Denom.Grado; Racional temp = Racional(ordsup, ordinf); temp.Numer = (a.Numer * b.Denom) + (a.Denom * b.Numer); temp.Denom = a.Denom * b.Denom; //Racional ress = Rectificar_Rac(temp); //return ress; return temp;
81
//Operador resta para estructuras de tipo Racional Racional operator -(Racional a,Racional b) int ordsup = 0; int ordinf = 0; if((a.Numer.Grado + b.Denom.Grado) >= (a.Denom.Grado + b.Numer.Grado)) ordsup = a.Numer.Grado + b.Denom.Grado; else ordsup = a.Denom.Grado + b.Numer.Grado; ordinf = a.Denom.Grado + b.Denom.Grado; Racional temp = Racional(ordsup, ordinf); temp.Numer = (a.Numer * b.Denom) - (a.Denom * b.Numer); temp.Denom = a.Denom * b.Denom; Racional ress = Rectificar_Rac(temp); return ress; //Operador Multiplicacion para estructuras de tipo Racional Racional operator *(Racional a,Racional b) int ordsup = 0; int ordinf = 0; ordsup = a.Numer.Grado + b.Numer.Grado; ordinf = a.Denom.Grado + b.Denom.Grado; Racional temp = Racional(ordsup, ordinf); temp.Numer = a.Numer * b.Numer; temp.Denom = a.Denom * b.Denom; Racional ress = Rectificar_Rac(temp); return ress; //Operador division para estructuras de tipo Racional Racional operator /(Racional a,Racional b) int ordsup = 0; int ordinf = 0; ordsup = a.Numer.Grado + b.Denom.Grado; ordinf = a.Denom.Grado + b.Numer.Grado; Racional temp = Racional(ordsup, ordinf); temp.Numer = a.Numer * b.Denom; temp.Denom = a.Denom * b.Numer; Racional ress = Rectificar_Rac(temp); return ress; //Fin operadores.
82
Anexo B. Código fuente de la variación para coeficientes reales del algoritmo de Jenkins-Truab para hallar las raíces de polinomios. /* rpoly.cpp -- Jenkins-Traub real polynomial root finder. * * (C) 2002, C. Bond. All rights reserved. * * Translation of TOMS493 from FORTRAN to C. This * implementation of Jenkins-Traub partially adapts * the original code to a C environment by restruction * many of the 'goto' controls to better fit a block * structured form. It also eliminates the global memory * allocation in favor of local, dynamic memory management. * * The calling conventions are slightly modified to return * the number of roots found as the function value. * * INPUT: * op - double precision vector of coefficients in order of * decreasing powers. * degree - integer degree of polynomial * * OUTPUT: * zeror,zeroi - output double precision vectors of the * real and imaginary parts of the zeros. * * RETURN: * returnval: -1 if leading coefficient is zero, otherwise * number of roots found. */ #include <stdio.h> #include <time.h> #include <math.h> #include <stdlib.h> void quad(double a,double b1,double c,double *sr,double *si, double *lr,double *li); void fxshfr(int l2, int *nz); void quadit(double *uu,double *vv,int *nz); void realit(double sss, int *nz, int *iflag); void calcsc(int *type); void nextk(int *type); void newest(int type,double *uu,double *vv); void quadsd(int n,double *u,double *v,double *p,double *q, double *a,double *b); double *p,*qp,*k,*qk,*svk;
83
double sr,si,u,v,a,b,c,d,a1,a2; double a3,a6,a7,e,f,g,h,szr,szi,lzr,lzi; double eta,are,mre; int n,nn,nmi,zerok; static int itercnt; int rpoly(double *op, int degree, double *zeror, double *zeroi, int info[] ) double t,aa,bb,cc,*temp,factor,rot; double *pt; double lo,max,min,xx,yy,cosr,sinr,xxx,x,sc,bnd; double xm,ff,df,dx,infin,smalno,base; int cnt,nz,i,j,jj,l,nm1,zerok; long sec; sec = clock(); /* The following statements set machine constants. */ base = 2.0; eta = 2.22e-16; infin = 3.4e38; smalno = 1.2e-38; are = eta; mre = eta; lo = smalno/eta; /* Initialization of constants for shift rotation. */ xx = sqrt(0.5); yy = -xx; rot = 94.0; rot *= 0.017453293; cosr = cos(rot); sinr = sin(rot); n = degree; /* Algorithm fails of the leading coefficient is zero. */ if (op[0] == 0.0) return -1; /* Remove the zeros at the origin, if any. */ while (op[n] == 0.0) j = degree - n; zeror[j] = 0.0; zeroi[j] = 0.0; n--; if (n < 1) return degree; /* * Allocate memory here */
84
temp = new double [degree+1]; pt = new double [degree+1]; p = new double [degree+1]; qp = new double [degree+1]; k = new double [degree+1]; qk = new double [degree+1]; svk = new double [degree+1]; /* Make a copy of the coefficients. */ for (i=0;i<=n;i++) p[i] = op[i]; /* Start the algorithm for one zero. */ _40: itercnt = 0; if (n == 1) zeror[degree-1] = -p[1]/p[0]; zeroi[degree-1] = 0.0; n -= 1; if( info != NULL ) info[ degree ] = 0; goto _99; /* Calculate the final zero or pair of zeros. */ if (n == 2) quad(p[0],p[1],p[2],&zeror[degree-2],&zeroi[degree-2], &zeror[degree-1],&zeroi[degree-1]); n -= 2; if( info != NULL ) info[ degree ] = info[ degree - 1] = 0; goto _99; /* Find largest and smallest moduli of coefficients. */ max = 0.0; min = infin; for (i=0;i<=n;i++) x = fabs(p[i]); if (x > max) max = x; if (x != 0.0 && x < min) min = x; /* Scale if there are large or very small coefficients. * Computes a scale factor to multiply the coefficients of the * polynomial. The scaling si done to avoid overflow and to * avoid undetected underflow interfering with the convergence * criterion. The factor is a power of the base. */ sc = lo/min; if (sc > 1.0 && infin/sc < max) goto _110;
85
if (sc <= 1.0) if (max < 10.0) goto _110; if (sc == 0.0) sc = smalno; l = (int)(log(sc)/log(base) + 0.5); factor = pow(base*1.0,l); if (factor != 1.0) for (i=0;i<=n;i++) p[i] = factor*p[i]; /* Scale polynomial. */ _110: /* Compute lower bound on moduli of roots. */ for (i=0;i<=n;i++) pt[i] = (fabs(p[i])); pt[n] = - pt[n]; /* Compute upper estimate of bound. */ x = exp((log(-pt[n])-log(pt[0])) / (double)n); /* If Newton step at the origin is better, use it. */ if (pt[n-1] != 0.0) xm = -pt[n]/pt[n-1]; if (xm < x) x = xm; /* Chop the interval (0,x) until ff <= 0 */ while (1) xm = x*0.1; ff = pt[0]; for (i=1;i<=n;i++) ff = ff*xm + pt[i]; if (ff <= 0.0) break; x = xm; dx = x; /* Do Newton interation until x converges to two * decimal places. */ while (fabs(dx/x) > 0.005) ff = pt[0]; df = ff; for (i=1;i<n;i++) ff = ff*x + pt[i]; df = df*x + ff; ff = ff*x + pt[n]; dx = ff/df; x -= dx;
86
itercnt++; bnd = x; /* Compute the derivative as the initial k polynomial * and do 5 steps with no shift. */ nm1 = n - 1; for (i=1;i<n;i++) k[i] = (double)(n-i)*p[i]/(double)n; k[0] = p[0]; aa = p[n]; bb = p[n-1]; zerok = (k[n-1] == 0); for(jj=0;jj<5;jj++) itercnt++; cc = k[n-1]; if (!zerok) /* Use a scaled form of recurrence if value of k at 0 is nonzero. */ t = -aa/cc; for (i=0;i<nm1;i++) j = n-i-1; k[j] = t*k[j-1]+p[j]; k[0] = p[0]; zerok = (fabs(k[n-1]) <= fabs(bb)*eta*10.0); else /* Use unscaled form of recurrence. */ for (i=0;i<nm1;i++) j = n-i-1; k[j] = k[j-1]; k[0] = 0.0; zerok = (k[n-1] == 0.0); /* Save k for restarts with new shifts. */ for (i=0;i<n;i++) temp[i] = k[i]; /* Loop to select the quadratic corresponding to each new shift. */ for (cnt = 0;cnt < 20;cnt++) /* Quadratic corresponds to a double shift to a * non-real point and its complex conjugate. The point * has modulus bnd and amplitude rotated by 94 degrees * from the previous shift. */ xxx = cosr*xx - sinr*yy;
87
yy = sinr*xx + cosr*yy; xx = xxx; sr = bnd*xx; si = bnd*yy; u = -2.0 * sr; v = bnd; fxshfr(20*(cnt+1),&nz); if (nz != 0) /* The second stage jumps directly to one of the third * stage iterations and returns here if successful. * Deflate the polynomial, store the zero or zeros and * return to the main algorithm. */ j = degree - n; zeror[j] = szr; zeroi[j] = szi; if( info != NULL ) info[ j + 1 ] = itercnt; n -= nz; for (i=0;i<=n;i++) p[i] = qp[i]; if (nz != 1) zeror[j+1] = lzr; zeroi[j+1] = lzi; if( info != NULL ) info[ j + 2 ] = 0; goto _40; /* If the iteration is unsuccessful another quadratic * is chosen after restoring k. */ for (i=0;i<n;i++) k[i] = temp[i]; /* Return with failure if no convergence after 20 shifts. */ _99: delete [] svk; delete [] qk; delete [] k; delete [] qp; delete [] p; delete [] pt; delete [] temp; info[ 0 ] = clock() - sec;
88
info[ 0 ] *= 1000; info[ 0 ] /= CLOCKS_PER_SEC; return degree - n; /* Computes up to L2 fixed shift k-polynomials, * testing for convergence in the linear or quadratic * case. Initiates one of the variable shift * iterations and returns with the number of zeros * found. */ void fxshfr(int l2,int *nz) double svu,svv,ui,vi,s; double betas,betav,oss,ovv,ss,vv,ts,tv; double ots,otv,tvv,tss; int type, i,j,iflag,vpass,spass,vtry,stry; *nz = 0; betav = 0.25; betas = 0.25; oss = sr; ovv = v; /* Evaluate polynomial by synthetic division. */ quadsd(n,&u,&v,p,qp,&a,&b); calcsc(&type); for (j=0;j<l2;j++) /* Calculate next k polynomial and estimate v. */ nextk(&type); calcsc(&type); newest(type,&ui,&vi); vv = vi; /* Estimate s. */ ss = 0.0; if (k[n-1] != 0.0) ss = -p[n]/k[n-1]; tv = 1.0; ts = 1.0; if (j == 0 || type == 3) goto _70; /* Compute relative measures of convergence of s and v sequences. */ if (vv != 0.0) tv = fabs((vv-ovv)/vv); if (ss != 0.0) ts = fabs((ss-oss)/ss); /* If decreasing, multiply two most recent convergence measures. */ tvv = 1.0; if (tv < otv) tvv = tv*otv; tss = 1.0; if (ts < ots) tss = ts*ots; /* Compare with convergence criteria. */
89
vpass = (tvv < betav); spass = (tss < betas); if (!(spass || vpass)) goto _70; /* At least one sequence has passed the convergence test. * Store variables before iterating. */ svu = u; svv = v; for (i=0;i<n;i++) svk[i] = k[i]; s = ss; /* Choose iteration according to the fastest converging * sequence. */ vtry = 0; stry = 0; if (spass && (!vpass) || tss < tvv) goto _40; _20: quadit(&ui,&vi,nz); if (*nz > 0) return; /* Quadratic iteration has failed. Flag that it has * been tried and decrease the convergence criterion. */ vtry = 1; betav *= 0.25; /* Try linear iteration if it has not been tried and * the S sequence is converging. */ if (stry || !spass) goto _50; for (i=0;i<n;i++) k[i] = svk[i]; _40: realit(s,nz,&iflag); if (*nz > 0) return; /* Linear iteration has failed. Flag that it has been * tried and decrease the convergence criterion. */ stry = 1; betas *=0.25; if (iflag == 0) goto _50; /* If linear iteration signals an almost double real * zero attempt quadratic iteration. */ ui = -(s+s); vi = s*s;
90
goto _20; /* Restore variables. */ _50: u = svu; v = svv; for (i=0;i<n;i++) k[i] = svk[i]; /* Try quadratic iteration if it has not been tried * and the V sequence is convergin. */ if (vpass && !vtry) goto _20; /* Recompute QP and scalar values to continue the * second stage. */ quadsd(n,&u,&v,p,qp,&a,&b); calcsc(&type); _70: ovv = vv; oss = ss; otv = tv; ots = ts; /* Variable-shift k-polynomial iteration for a * quadratic factor converges only if the zeros are * equimodular or nearly so. * uu, vv - coefficients of starting quadratic. * nz - number of zeros found. */ void quadit(double *uu,double *vv,int *nz) double ui,vi; double mp,omp,ee,relstp,t,zm; int type,i,j,tried; *nz = 0; tried = 0; u = *uu; v = *vv; j = 0; /* Main loop. */ _10: itercnt++; quad(1.0,u,v,&szr,&szi,&lzr,&lzi); /* Return if roots of the quadratic are real and not * close to multiple or nearly equal and of opposite
91
* sign. */ if (fabs(fabs(szr)-fabs(lzr)) > 0.01 * fabs(lzr)) return; /* Evaluate polynomial by quadratic synthetic division. */ quadsd(n,&u,&v,p,qp,&a,&b); mp = fabs(a-szr*b) + fabs(szi*b); /* Compute a rigorous bound on the rounding error in * evaluating p. */ zm = sqrt(fabs(v)); ee = 2.0*fabs(qp[0]); t = -szr*b; for (i=1;i<n;i++) ee = ee*zm + fabs(qp[i]); ee = ee*zm + fabs(a+t); ee *= (5.0 *mre + 4.0*are); ee = ee - (5.0*mre+2.0*are)*(fabs(a+t)+fabs(b)*zm); ee = ee + 2.0*are*fabs(t); /* Iteration has converged sufficiently if the * polynomial value is less than 20 times this bound. */ if (mp <= 20.0*ee) *nz = 2; return; j++; /* Stop iteration after 20 steps. */ if (j > 20) return; if (j < 2) goto _50; if (relstp > 0.01 || mp < omp || tried) goto _50; /* A cluster appears to be stalling the convergence. * Five fixed shift steps are taken with a u,v close * to the cluster. */ if (relstp < eta) relstp = eta; relstp = sqrt(relstp); u = u - u*relstp; v = v + v*relstp; quadsd(n,&u,&v,p,qp,&a,&b); for (i=0;i<5;i++) calcsc(&type); nextk(&type); tried = 1; j = 0; _50:
92
omp = mp; /* Calculate next k polynomial and new u and v. */ calcsc(&type); nextk(&type); calcsc(&type); newest(type,&ui,&vi); /* If vi is zero the iteration is not converging. */ if (vi == 0.0) return; relstp = fabs((vi-v)/vi); u = ui; v = vi; goto _10; /* Variable-shift H polynomial iteration for a real zero. * sss - starting iterate * nz - number of zeros found * iflag - flag to indicate a pair of zeros near real axis. */ void realit(double sss, int *nz, int *iflag) double pv,kv,t,s; double ms,mp,omp,ee; int i,j; *nz = 0; s = sss; *iflag = 0; j = 0; /* Main loop */ while (1) itercnt++; pv = p[0]; /* Evaluate p at s. */ qp[0] = pv; for (i=1;i<=n;i++) pv = pv*s + p[i]; qp[i] = pv; mp = fabs(pv); /* Compute a rigorous bound on the error in evaluating p. */ ms = fabs(s); ee = (mre/(are+mre))*fabs(qp[0]); for (i=1;i<=n;i++) ee = ee*ms + fabs(qp[i]); /* Iteration has converged sufficiently if the polynomial * value is less than 20 times this bound.
93
*/ if (mp <= 20.0*((are+mre)*ee-mre*mp)) *nz = 1; szr = s; szi = 0.0; return ; // HVE return added j++; /* Stop iteration after 10 steps. */ if (j > 10) return; if (j < 2) goto _50; if (fabs(t) > 0.001*fabs(s-t) || mp < omp) goto _50; /* A cluster of zeros near the real axis has been * encountered. Return with iflag set to initiate a * quadratic iteration. */ *iflag = 1; sss =s; // HVE sss=s added return; /* Return if the polynomial value has increased significantly. */ _50: omp = mp; /* Compute t, the next polynomial, and the new iterate. */ kv = k[0]; qk[0] = kv; for (i=1;i<n;i++) kv = kv*s + k[i]; qk[i] = kv; if (fabs(kv) <= fabs(k[n-1])*10.0*eta) // HVE n -> n-1 /* Use unscaled form. */ k[0] = 0.0; for (i=1;i<n;i++) k[i] = qk[i-1]; else /* Use the scaled form of the recurrence if the value * of k at s is nonzero. */ t = -pv/kv; k[0] = qp[0]; for (i=1;i<n;i++) k[i] = t*qk[i-1] + qp[i]; kv = k[0]; for (i=1;i<n;i++) kv = kv*s + k[i];
94
t = 0.0; if (fabs(kv) > fabs(k[n-1]*10.0*eta)) t = -pv/kv; s += t; /* This routine calculates scalar quantities used to * compute the next k polynomial and new estimates of * the quadratic coefficients. * type - integer variable set here indicating how the * calculations are normalized to avoid overflow. */ void calcsc(int *type) /* Synthetic division of k by the quadratic 1,u,v */ quadsd(n-1,&u,&v,k,qk,&c,&d); if (fabs(c) > fabs(k[n-1]*100.0*eta)) goto _10; if (fabs(d) > fabs(k[n-2]*100.0*eta)) goto _10; *type = 3; /* Type=3 indicates the quadratic is almost a factor of k. */ return; _10: if (fabs(d) < fabs(c)) *type = 1; /* Type=1 indicates that all formulas are divided by c. */ e = a/c; f = d/c; g = u*e; h = v*b; a3 = a*e + (h/c+g)*b; a1 = b - a*(d/c); a7 = a + g*d + h*f; return; *type = 2; /* Type=2 indicates that all formulas are divided by d. */ e = a/d; f = c/d; g = u*b; h = v*b; a3 = (a+g)*e + h*(b/d); a1 = b*f-a; a7 = (f+u)*a + h; /* Computes the next k polynomials using scalars * computed in calcsc.
95
*/ void nextk(int *type) double temp; int i; if (*type == 3) /* Use unscaled form of the recurrence if type is 3. */ k[0] = 0.0; k[1] = 0.0; for (i=2;i<n;i++) k[i] = qk[i-2]; return; temp = a; if (*type == 1) temp = b; if (fabs(a1) <= fabs(temp)*eta*10.0) /* If a1 is nearly zero then use a special form of the * recurrence. */ k[0] = 0.0; k[1] = -a7*qp[0]; for(i=2;i<n;i++) k[i] = a3*qk[i-2] - a7*qp[i-1]; return; // HVE return added /* Use scaled form of the recurrence. */ a7 /= a1; a3 /= a1; k[0] = qp[0]; k[1] = qp[1] - a7*qp[0]; for (i=2;i<n;i++) k[i] = a3*qk[i-2] - a7*qp[i-1] + qp[i]; /* Compute new estimates of the quadratic coefficients * using the scalars computed in calcsc. */ void newest(int type,double *uu,double *vv) double a4,a5,b1,b2,c1,c2,c3,c4,temp; /* Use formulas appropriate to setting of type. */ if (type == 3) /* If type=3 the quadratic is zeroed. */ *uu = 0.0;
96
*vv = 0.0; return; if (type == 2) a4 = (a+g)*f + h; a5 = (f+u)*c + v*d; else a4 = a + u*b +h*f; a5 = c + (u+v*f)*d; /* Evaluate new quadratic coefficients. */ b1 = -k[n-1]/p[n]; b2 = -(k[n-2]+b1*p[n-1])/p[n]; c1 = v*b2*a1; c2 = b1*a7; c3 = b1*b1*a3; c4 = c1 - c2 - c3; temp = a5 + b1*a4 - c4; if (temp == 0.0) *uu = 0.0; *vv = 0.0; return; *uu = u - (u*(c3+c2)+v*(b1*a1+b2*a7))/temp; *vv = v*(1.0+c4/temp); return; /* Divides p by the quadratic 1,u,v placing the quotient * in q and the remainder in a,b. */ void quadsd(int nn,double *u,double *v,double *p,double *q, double *a,double *b) double c; int i; *b = p[0]; q[0] = *b; *a = p[1] - (*b)*(*u); q[1] = *a; for (i=2;i<=nn;i++) c = p[i] - (*a)*(*u) - (*b)*(*v); q[i] = c; *b = *a; *a = c;
97
/* Calculate the zeros of the quadratic a*z^2 + b1*z + c. * The quadratic formula, modified to avoid overflow, is used * to find the larger zero if the zeros are real and both * are complex. The smaller real zero is found directly from * the product of the zeros c/a. */ void quad(double a,double b1,double c,double *sr,double *si, double *lr,double *li) double b,d,e; if (a == 0.0) /* less than two roots */ if (b1 != 0.0) *sr = -c/b1; else *sr = 0.0; *lr = 0.0; *si = 0.0; *li = 0.0; return; if (c == 0.0) /* one real root, one zero root */ *sr = 0.0; *lr = -b1/a; *si = 0.0; *li = 0.0; return; /* Compute discriminant avoiding overflow. */ b = b1/2.0; if (fabs(b) < fabs(c)) if (c < 0.0) e = -a; else e = a; e = b*(b/fabs(c)) - e; d = sqrt(fabs(e))*sqrt(fabs(c)); else e = 1.0 - (a/b)*(c/b); d = sqrt(fabs(e))*fabs(b); if (e < 0.0) /* complex conjugate zeros */ *sr = -b/a; *lr = *sr; *si = fabs(d/a); *li = -(*si); else
98
if (b >= 0.0) d = -d; /* real zeros. */ *lr = (-b+d)/a; *sr = 0.0; if (*lr != 0.0) *sr = (c/ *lr)/a; *si = 0.0; *li = 0.0;
99
Anexo C. Código fuente de la configuración del hardware de la ADSP BF533 en VisualDSP++. El siguiente código indica el orden en que las funciones deben ser llamadas:
Init_EBIU(); Init_LED(); Init_Flash(); Init1836(); Init_Sport0(); Init_DMA(); Init_Interrupts(); Enable_DMA_Sport0(); Abajo el código fuente de cada una de las funciones utilizadas arriba: //--------------------------------------------------------------------------// // Funcion: Init_EBIU // // Descripcion: Esta funcion inicializa y habilita los bancos de memoria // asíncrona en la EBIU para que la Flash A pueda ser accesada //--------------------------------------------------------------------------// void Init_EBIU(void) *pEBIU_AMBCTL0 = 0x7bb07bb0; *pEBIU_AMBCTL1 = 0x7bb07bb0; *pEBIU_AMGCTL = 0x000f; //--------------------------------------------------------------------------// // Funcion: Init_Flash // // Descripcion: Esta funcion inicializa la direccion de Port A en la Flash A // Se ubica en el ultimo bit del control Flash A port A para // resetearl //--------------------------------------------------------------------------// void Init_Flash(void) *pFlashA_PortA_Dir = 0x1; //--------------------------------------------------------------------------// // Funcion: Init1836 // // Descripcion: Esta funcion configura el codec AD1836 // El contenido del arreglo Codec1836TxRegs es enviado al
100
// codec. //--------------------------------------------------------------------------// void Init1836(void) int i; int j; static unsigned char ucActive_LED = 0x01; *pFlashA_PortA_Data = 0x00; *pFlashA_PortA_Data = ucActive_LED; for (i=0; i<0xf0000; i++) asm("nop;"); *pSPI_FLG = FLS4; *pSPI_BAUD = 16; *pSPI_CTL = TIMOD_DMA_TX | SIZE | MSTR; *pDMA5_PERIPHERAL_MAP = 0x5000; *pDMA5_CONFIG = WDSIZE_16; *pDMA5_START_ADDR = (void *)sCodec1836TxRegs; *pDMA5_X_COUNT = CODEC_1836_REGS_LENGTH; *pDMA5_X_MODIFY = 2; *pDMA5_CONFIG = (*pDMA5_CONFIG | DMAEN); *pSPI_CTL = (*pSPI_CTL | SPE); for (j=0; j<0xaff0; j++) asm("nop;"); *pSPI_CTL = 0x0000; //--------------------------------------------------------------------------// // Funcion: Init_Sport0 // // Descripcion: Configurar el puerto Sport0 para el modo I2S, recibir y // transmitir datos desde y hacia el AD1836. // Configurar el Sport para relog externo y sincronizar tramas. //--------------------------------------------------------------------------// void Init_Sport0(void) *pSPORT0_RCR1 = RFSR | RCKFE; *pSPORT0_RCR2 = SLEN_24 | RXSE | RSFSE;
101
*pSPORT0_TCR1 = TFSR | TCKFE; *pSPORT0_TCR2 = SLEN_24 | TXSE | TSFSE; //-------------------------------------------------------------------------- // // Funcion: Init_DMA // // Descripcion: Configurar la interrupcion del DMA1 y DMA2 para el control // de registros recepcion y trasmicion desde y hacia el SPORT0 // El canal 1 Y 2 del DMA esta mapeado por defecto al SPORT0 // RX y SPORTO TX respectivamente //--------------------------------------------------------------------------// void Init_DMA(void) *pDMA1_PERIPHERAL_MAP = 0x1000; // Recordemos que el tamaño de bits a procesar es de 32, por tanto WDSIZE=32 *pDMA1_CONFIG = WNR | WDSIZE_32 | DI_EN | FLOW_1; *pDMA1_START_ADDR = (void *)iRxBuffer1; *pDMA1_X_COUNT = 2; *pDMA1_X_MODIFY = 4; *pDMA2_PERIPHERAL_MAP = 0x2000; *pDMA2_CONFIG = WDSIZE_32 | FLOW_1; *pDMA2_START_ADDR = (void *)iTxBuffer1; *pDMA2_X_COUNT = 2; *pDMA2_X_MODIFY = 4; //--------------------------------------------------------------------------// // Funcion: Enable_DMA_Sport // // Descripcion: Habilitar DMA1, DMA2, Sport0 TX y Sport0 RX //--------------------------------------------------------------------------// void Enable_DMA_Sport0(void) *pDMA2_CONFIG = (*pDMA2_CONFIG | DMAEN); *pDMA1_CONFIG = (*pDMA1_CONFIG | DMAEN); *pSPORT0_TCR1 = (*pSPORT0_TCR1 | TSPEN); *pSPORT0_RCR1 = (*pSPORT0_RCR1 | RSPEN);
102
//--------------------------------------------------------------------------// // Function: Init_Interrupts // // Description: Inicializar interrupciones para Sport0 RX //--------------------------------------------------------------------------// void Init_Interrupts(void) *pSIC_IAR0 = 0xffffffff; *pSIC_IAR1 = 0xffffff2f; *pSIC_IAR2 = 0xffffffff; register_handler(ik_ivg9, Sport0_RX_ISR); *pSIC_IMASK = 0x00000200; //--------------------------------------------------------------------------// // Function: Init_LED // // Description: Inicializa el uso de los LED del DSP //--------------------------------------------------------------------------// void Init_LED(void) *pFlashA_PortB_Dir = 0x3f; *pFlashA_PortB_Data = 0x3f; //Interrupt que sucede cuando el codec entrega una nueva muestra EX_INTERRUPT_HANDLER(Sport0_RX_ISR) // Se confirma la interrupción *pDMA1_IRQ_STATUS = 0x0001; iChannel0LeftIn = iRxBuffer1[INTERNAL_ADC_L0]; iChannel0RightIn = iRxBuffer1[INTERNAL_ADC_R0]; iChannel1LeftIn = iRxBuffer1[INTERNAL_ADC_L1]; iChannel1RightIn = iRxBuffer1[INTERNAL_ADC_R1];
Process_Data(); //Filtrado de la señal iTxBuffer1[INTERNAL_DAC_L0] = iChannel0LeftOut; iTxBuffer1[INTERNAL_DAC_R0] = iChannel0RightOut; iTxBuffer1[INTERNAL_DAC_L1] = iChannel1LeftOut; iTxBuffer1[INTERNAL_DAC_R1] = iChannel1RightOut;
