Las investigaciones realizadas en el campo de la medicina dieron lugar en 1977 a la concesión del Premio Nobel a Sir Godfrey N. y Allan M. Cormack, por la invención del sistema de exploración mediante tomografía axial computerizada (T.A.C.), revolucionando a partir de ese momento, el diagnóstico en el campo de la radiología.
Actualmente la tomográfica computerizada es una práctica rutinaria en la mayoría de los hospitales y la investigación está centrada en tres aspectos:
La necesidad de obtener resoluciones con una mayor calidad tiene dos implicaciones directas, el aumento en la resolución con la que se trabaje (mayor nº de puntos a tratar) y como consecuencia el aumento del tiempo de computación. La utilización por tanto de supercomputadores está sobradamente justificada en este campo dada la cantidad de datos a tratar y el coste computacional de los algoritmos de reconstrucción, como veremos posteriormente.
En los últimos años, la aparición de nuevas técnicas en las que se utiliza radiación-gamma, conocida generalmente como tomografía de Emisión (PET y SPECT), ha permitido conocer, no solo la estructura interna de los cuerpos, sino el aspecto funcional de determinados órganos o tejidos.
Centraremos en este artículo nuestra atención en la aceleración del tiempo de reconstrucción, exponiendo en primer lugar los fundamentos científicos que rigen esta técnica, para después analizar los problemas que se presentan al abordar su tratamiento desde el punto de vista computacional y los sistemas a utilizar.
¿QUE RECONSTRUYO Y A PARTIR DE QUE DATOS?
Si iluminamos un objeto con un haz de rayos-X es sobradamente conocido que obtenemos una imagen bidimensional, la cual nos ensaña la silueta de la proyección objeto sobre el plano de recepción de la radiación junto a un conjunto de características interiores del objeto, cuya visualización depende de la mayor o menor absorción que determinadas zonas del objeto tengan al paso de los rayos-X.
Cuando reconstruimos estamos realizando el proceso inverso al descrito anteriormente, es decir, que a partir de un conjunto de estas proyecciones, obtenidas mediante determinadas geometrías de proyección, conseguimos reconstruir el objeto que dio lugar a tales proyecciones.
Modelizando idealmente el problema, la ecuación que rige el fenómeno de proyectar un objeto es la siguiente:
N(x) = No e-µx donde No es el nº de fotones incidentes en un intervalo de tiempo, N(x) el nº de fotones a una distancia x del borde del objeto y µ es el coeficiente de atenuación del material (considerado constante en esta primera aproximación).
Si consideramos a µ dependiente de dos coordenadas espaciales (x, y), la ecuación que nos indicaría el nº de fotones salientes de un objeto (Nd) en relación al nº de fotones incidentes para un solo rayo vendrá expresado entonces de la siguiente forma:
o la equivalente:
El lado izquierdo de esta ecuación constituye la integral de un rayo de una proyección. De esta forma medidas de (Nin/Nd) tomadas para diferentes rayos y diferentes ángulos (geometrías determinadas), pueden ser interpretadas como un conjunto de proyecciones de la función µ(x,y).
Podemos concluir, por tanto, que nuestro conjunto de datos lo constituye un conjunto mínimo de proyecciones y el objeto a reconstruir vendrá identificado por la función µ(x,y,z) para 3D.
Un ejemplo de reconstrucción es el que se observa en las figuras1 y 2, en las cuales se representa respectivamente el modelo conocido como fantasma de Logan (ampliamente utilizado en tomografía para comparar la calidad de reconstrucción de distintos algoritmos) y una reconstrucción real.
DIFERENTES ALGORITMOS
La solución algorítmica a este problema sigue dos caminos diferentes, uno de ellos es el de los métodos directos o transformados y el otro el de los métodos algebraicos de reconstrucción. El más significativo de los primeros, ampliamente utilizado en la mayoría de los tomografos actuales, es al denominado Filtered Back Projection (FPB), basado en el teorema de Fourier de la proyección central y tiene los siguientes pasos.
Para estos métodos el problema se plantea como el de encontrar la solución al siguiente sistema de ecuaciones lineales:
g = H f + n donde g y n representan los vectores de datos y ruido respectivamente, f el vector de incógnitas (µ) y H es la matriz de transformación asociada al modelo físico del proceso de adquisición de datos. Cada elemento Hij de la matriz establece una relación entre el rayo i y el elemento de volumen j (voxel).
El algoritmo emblemático de esta familia es el conocido ART, el cual actualiza de forma interactiva cada una de las incógnitas del sistema planteado utilizando la siguiente ecuación.
En donde k es el índice de iteraciones.
Siendo N la dimensión del problema (volumen a reconstruir N3 y proyecciones N2) y M el nº de proyecciones utilizados, tendremos una complejidad algorítmica del orden O(MN5). Si consideramos para casos reales N=128 y M=300 nos encontramos con un problema cuyo coste computacional es muy elevado. (221 incógnitas y 300 x 214 ecuaciones) al cual hay que añadir la también elevada necesidad de memoria. (8 Mbyte solo para las incógnitas).
El dar una medida de tiempo en el cual obtengamos un resultado, como es lógico dependerá de la máquina utilizada, pero en cualquier caso estamos hablando de horas para el caso de máquinas monoprocesador, cuyo ejemplo podemos indicar los 241 minutos empleados por una Sparstation 2, para una sola iteración (la reconstrucción requiere bastantes más iteraciones hasta alcanzar una convergencia adecuada).
Todo esto nos indica que una solución eficaz a este problema requiere la utilización de máquinas que tengan una elevada capacidad de cálculo, justificándose el empleo de computadores masivamente paralelos.
PARALELISMO Y SUPERCOMPUTADORES
El término supercomputador se asocia a aquellos computadores que han sido diseñados para poder ejecutar un gran número de operaciones por segundo, sobre números codificados normalmente en coma flotante, y en cada época, son los que efectúan un mayor número de operaciones por segundo.
A principios de los 70, la velocidad de estas máquinas se situaba por debajo de los 100 MFlopx (108 operaciones/seg). En la actualidad, existen supercomputadores con velocidades de pico por encima de los 100 GFlops (1011 operaciones/seg). Esta velocidad junto a la enorme capacidad de las memorias y de la mejora en los sistemas de entrada/salida, ha permitido abordar los considerados grandes retos, como son: modelado del clima, turbulencia de fluidos y diseño de nuevos materiales.
Está bastante admitido que las arquitecturas futuras para supercomputación estarán basadas en sistemas multiprocesadores (Intel, Kendall, Square, Alliant, Alpha,...) con memoria local DRAM conectados a una Red de Interconexión y con interfaces de red entre los pares memoria-procesador que permitan al usuario adaptar la topología al programa particular que está corriendo. Es esta tecnología fue pionera Meiko en sus máquinas Computer Surface basadas en trasputers estándar de cuatro links, utilizando un interface de programación denominado CHIMP. Ejemplos actuales de esta organización, esquematicamente representada en la figura 2, los tenemos en el Cray Research T3D, Intel Delta, Parangon y CM-5. La analogía de los nodos de procesamiento paralelo y estaciones de trabajo sugiere que el procesamiento paralelo del futuro pueda utilizar computadores comerciales (CPUs de propósito general)
Si las redes de área local de alta velocidad conectan computadores de mesa a través de un conmutador central de gran anchura de banda y los gastos de envío del mensaje pueden reducirse en las estaciones de trabajo, entonces la distribución entre procesado-res paralelos y «clusters» de estaciones de trabajo puede desaparecer.
Esta breve introducción a los supercompu-tadores nos sirve para justificar la utilización de sistemas multipro-cesador MIMD con memoria distribuida para resolver el sistema de ecuaciones representativo del problema de reconstrucción de imágenes.
Como es sabido la penalización que sufrimos utilizando estos sistemas es la necesaria intercomunicación de datos entre procesadores y como consecuencia el descenso del rendimiento obtenido. Con el fin de evitar tal descenso de forma significativa debemos resolver los problemas típicos del procesamiento paralelo:
El modelo de programación SPMD dado a los algoritmos de tratamiento de imágenes, caracterizados por operaciones locales sobre vecinos próximos, puede aplicarse también al problema de reconstrucción obteniendo resultados con rendimientos superiores al 80%, dependiendo del tamaño del problema (resolución en la que se trabaje) y del número de procesadores utilizado (el rendimiento puede descender si utilizamos un número excesivo de procesadores).
La utilización de nuevas geometrías de proyección como cone-beam, junto con la aparición de sistemas en los que la potencia de cálculo por nodo se conjugue de forma eficaz con la obtención de un alto rendimiento mediante su programación, hará que en un futuro próximo la calidad de reconstrucciones aumente de forma espectacular, facilitando con ello un mejor diagnóstico en medicina.