Desarrollo de una interfaz gráfica de usuario en entorno MATLAB para el manejo de software de segmentación de imágenes RMN

UNIVERSIDAD POLITÉCNICA DE VALENCIA ESCUELA TÉCNICA SUPERIOR DE INGENIERÍA INFORMÁTICA PROYECTO FIN DE CARRERA Desarrollo de una interfaz gráfica de usuario en entorno MATLAB para el manejo de software de segmentación de imágenes RMN Realizado por: Javier Cancio del Busto Dirigido por: José Vicente Manjón Herrera Valencia, Septiembre de 2010 “Investigar es ver lo que todo el mundo ha visto, y pensar lo que nadie más ha pensado." Albert Szent-Györgyi (Filósofo húngaro, 1893-1986) Índice de Contenidos 1. 1.1. 1. INTRODUCCIÓN INTRODUCCIÓNINTRODUCCIÓN INTRODUCCIÓN 8 88 8 1.1 P RÓLOGO 8 1.2 S ITUACIÓN Y OBJETIVOS 9 1.3 E STRUCTURA DEL DOCUMENTO 10 2. 2.2. 2. CONCEPTOS BÁSICOS DE CONCEPTOS BÁSICOS DECONCEPTOS BÁSICOS DE CONCEPTOS BÁSICOS DE RESONAN RESONAN RESONAN RESONANCIA MAGNÉTICA NUCLEAR 1 CIA MAGNÉTICA NUCLEAR 1CIA MAGNÉTICA NUCLEAR 1 CIA MAGNÉTICA NUCLEAR 11 11 1 2.1 ORIGEN Y APLICACIONES 11 2.2 PRICIPIOS FÍSICOS DE LA RMN 12 2.3 IMÁGENES MÉDICAS POR RMN 14 2.4 EL CONCEPTO DE VÓXEL 19 3. 3.3. 3. MATERIAL Y MÉTODOS MATERIAL Y MÉTODOSMATERIAL Y MÉTODOS MATERIAL Y MÉTODOS UTILIZADOS UTILIZADOS UTILIZADOS UTILIZADOS 20 20 20 20 3.1 IMÁGENES MÉDICAS 20 3.2 MATLAB 21 3.3 SPM 23 3.4 MÓDULOS 25 4. 4.4. 4. BRAINSEG BRAINSEGBRAINSEG BRAINSEG 28 28 28 28 5. 5.5. 5. RESULTADOS RESULTADOSRESULTADOS RESULTADOS 30 30 30 30 6. 6. 6. 6. DISC DISCDISC DISCUSIÓN Y CONCLUSIONES USIÓN Y CONCLUSIONES USIÓN Y CONCLUSIONES USIÓN Y CONCLUSIONES 33 33 33 33 7. 7.7. 7. AGRADECIMIENTOS AGRADECIMIENTOSAGRADECIMIENTOS AGRADECIMIENTOS 34 34 34 34 8 88 8. .. . BIBLIOGRAFÍA BIBLIOGRAFÍABIBLIOGRAFÍA BIBLIOGRAFÍA 3 3 3 35 55 5 ANEXO I ANEXO IANEXO I ANEXO I 36 36 36 36 Introducción 8 1. 1.1. 1. Introduc IntroducIntroduc Introducción ciónción ción 1.1 1.11.1 1.1 Prólogo PrólogoPrólogo Prólogo El inicio del siglo XIX vino marcado por el gran avance que supuso la “Revolución Industrial”. De forma análoga, durante el siglo XX se produce otro proceso revolucionario, gracias al avance de la ciencia, y concretamente de la “Ciencia de la computación”, a la que aún le queda un largo camino. Hoy en día, nos resulta prácticamente imposible pensar en la vida sin la ayuda que nos proporciona la informática en casi todos los ámbitos : económicos, sociales y tecnológicos. Uno de los campos afectados por este avance es el campo de la salud. Gracias al progreso tecnológico ha sido posible juntar disciplinas aparentemente divergentes como el de la informática y la medicina, y gracias a esta convergencia ha sido posible realizar diagnósticos con una velocidad mucho mayor, explorar las causas de un gran número de enfermedades y aumentar la eficacia de los tratamientos. El presente texto se centra en la integración de diferentes métodos tecnológicos, desarrollados para poder tratar con una mayor eficiencia las enfermedades que derivan en un cambio en el volumen de ciertos tejidos cerebrales, como es el caso del Alzheimer. El Alzheimer es una enfermedad neurodegenerativa que se muestra como un deterioro cognitivo y trastornos conductuales. Se caracteriza típicamente por una pérdida progresiva de memoria y otras capacidades mentales, a medida que las neuronas mueren y se atrofian diferentes zonas del cerebro. Es precisamente la muerte de las células nerviosas lo que causa pérdida de tejido, especialmente en la región del hipocampo, hasta en un 32% 1 . Es por estas enfermedades el interés de crear una herramienta gráfica que simplifique el uso de los métodos de segmentación en imágenes de resonancia magnética nuclear, junto con métodos previos, como el filtrado de imágenes y también poder procesar simultáneamente múltiples casos. 1 Olivier Colliot, Gaël Chételat, Marie Chupin. “Discrimination between Alzheimer Disease, Mild cognitive Impairment, and Normal Aging by Using Automated Segmentation of the Hippocampus” Introducción 9 1.2 1.21.2 1.2 Situación y objetivos Situación y objetivosSituación y objetivos Situación y objetivos El estudio actual se enmarca dentro del área de investigación del grupo ITACA (Instituto de Aplicaciones de las Tecnologías de la Información y de las Comunicaciones Avanzadas). El grupo ITACA cuenta con un total de cinco secciones de I+D+i entre las que se encuentra el grupo de informática biomédica (IBIME). Dentro de esta área existe un departamento que se ocupa del tratamiento informatizado de las imágenes médicas digitales facilitando las tareas de análisis, procesado, reconstrucción, segmentación, etc. Figura 1.1. Figura 1.1.Figura 1.1. Figura 1.1. Jerarquía del grupo de investigación IBIME. Las imágenes que se han utilizado en este proyecto corresponden al campo de la resonancia magnética nuclear (en adelante RM). Las continuas mejoras que se han producido en este tipo de imágenes durante los últimos años, han propiciado la ejecución de numerosos ensayos debido a la mayor resolución y nivel de precisión que aportan a los equipos. Dentro del abanico de posibilidades que ofrece la tecnología actual, la herramienta desarrollada trabajará con imágenes RM pertenecientes a una zona anatómica muy concreta: el encéfalo. El encéfalo se encuentra situado en el interior del cráneo y está constituido por el cerebro, el cerebelo y el tronco cerebral. El objetivo del presente proyecto es el desarrollo de una interfaz gráfica de usuario (GUI) que sirva de herramienta a los operarios de las imágenes RM para poder realizar las tareas de segmentación de forma más fácil y cómoda. La herramienta no sólo tiene como finalidad la segmentación de las Conceptos básicos de Resonancia Magnética Nuclear 16 • El paciente emite una señal que es recibida y utilizada para reconstruir la imagen ¿Por qué se coloca al paciente dentro de un campo magnético? Las partículas subatómicas, como hemos comentario en el punto anterior, tienen propiedades inherentes. Protones y electrones tienen carga mientras que los neutrones no. Del mismo modo que la carga de un átomo es la suma de todas sus cargas, el espín total será la suma de los espines de todas sus partículas subatómicas. Imaginemos el átomo más sencillo: un protón. La carga del protón realiza un movimiento de espín con él ¿y qué es una carga eléctrica en movimiento? Es una corriente eléctrica. Y las cargas eléctricas producen campos magnéticos. Así pues, ese protón (y como él cualquier núcleo) cargado y girando produce a su alrededor un campo magnético. Figura 2.4. Figura 2.4.Figura 2.4. Figura 2.4. Campo magnético creado por un protón. En la mayor parte de los materiales, los momentos magnéticos están orientados al azar (figura 2.2), anulándose unos a otros y dando como resultado un momento magnético total de cero. Pero, ¿qué ocurre cuando colocamos los átomos bajo un intenso campo magnético externo? Los momentos magnéticos se alinean en la dirección del campo magnético como pequeñas brújulas. Cierto que se alinean en la misma dirección pero pueden estar en sentidos contrarios. Aunque, en realidad, hay más protones que se alinean en el sentido de las líneas de campo magnético externo porque este estado tiene menos energía (de la misma manera que nosotros podemos alinearnos con las líneas del campo gravitatorio andando sobre nuestras manos o nuestros pies, este último estado es más estable energéticamente). Queda así un remanente de protones que no anulan sus momentos magnéticos con sus opuestos. En estas condiciones el paciente se transforma en un imán. Conceptos básicos de Resonancia Magnética Nuclear 17 Pero los protones no giran con un movimiento limpio alrededor de un eje perfectamente alineado con el campo magnético externo, sino que al igual que una peonza sometida a un campo gravitatorio describen un giro de forma cónica. Los protones sometidos a un campo magnético describen el mismo giro, llamado precesión . La rapidez con la que «precesa» un núcleo se llama frecuencia de precesión y depende exclusivamente de la fuerza del campo magnético externo y del tipo de núcleo implicado. Esta relación queda definida por la ecuación de Larmor : W = γ 0 B r Ecuaci EcuaciEcuaci Ecuación 2.2 ón 2.2ón 2.2 ón 2.2. .. . Ecuación de Larmor. En esta ecuación, W representa la frecuencia de precesión (en Hz), 0 B r la fuerza del campo magnético (en teslas) y γ la constante giromagnética (utilizada en la ecuación 2.1). Estas propiedades del campo magnético pueden utilizarse para analizar muestras. Supongamos que colocamos una muestra de material desconocido en un campo magnético de valor conocido: Si obtenemos un valor de la frecuencia de precesión podemos calcular la constante giromagnética, dado que es específica para cada material, sabremos de qué material se trata. La fuerza de la señal nos indicará la cantidad de material de ese tipo que hay en la muestra. Sin embargo, los núcleos precesan al azar, es decir, en un instante dado, existen núcleos orientados en cualquier dirección dentro de su giro de precesión. Fuerzas magnéticas orientadas en direcciones opuestas se cancelan unas a otras. Esto es así para todas las direcciones excepto para el eje Z (a lo largo del campo magnético externo). En la dirección Z los vectores suman sus fuerzas magnéticas dando como resultado final un vector magnético en la misma dirección que el campo externo. Figura 2.5. Figura 2.5.Figura 2.5. Figura 2.5. Momento magnético resultante de todos los µ r . Conceptos básicos de Resonancia Magnética Nuclear 18 En estas condiciones el paciente se transforma en un imán y este vector que emite se llama magnetización longitudinal . Es este nuevo vector el que puede ser utilizado para obtener la señal, el problema es que mientras esté paralelo a las líneas del campo magnético externo no podemos cuantificarlo. Parece que entonces no hemos conseguido nada, pero continuemos con el segundo paso: ondas de radio. Para campos magnéticos de intensidad razonablemente elevada, la frecuencia de precesión de casi todos los núcleos de interés queda dentro de la banda de radiofrecuencia (RF). Si mandamos a un paciente «imantado» un pulso de radiofrecuencia (PRF) igual a la frecuencia de precesión del hidrógeno bajo ese campo magnético se producirá un fenómeno conocido como resonancia. Se puede ilustrar con el uso de diapasones; si usted entra en una habitación con diapasones afinados en do, re, mi y toca un diapasón de frecuencia do, hará vibrar a todos los diapasones que estén en do y sólo a esos. Pues bien, los H + absorberán energía de esa radiofrecuencia y se pasarán a su estado de alta energía, es decir se colocarán en sentido opuesto al campo magnético externo, reduciéndose la magnetización longitudinal (recordamos que los momentos en sentidos opuestos se anulan). Pero, no sólo ocurre eso, sino que además los H + comenzarán a precesar en fase, por lo tanto ahora sus momentos magnéticos en las direcciones X e Y se suman en vez de anularse. Aparece entonces un vector de magnetización transversal. Este nuevo vector de magnetización transversal, al ser la suma de los momentos de cada protón, precesa de la misma forma. Si miramos desde fuera se acerca nosotros para luego alejarse y después volver a acercarse. Este vector magnético en movimiento induce una corriente eléctrica que puede ser registrada por una antena. Cuanto mayor sea el vector, mejor será la imagen de RM obtenida. El valor del mismo depende de γ, 0 B r y la densidad de espines. Bien, ya tenemos una señal, pero ¿de qué parte del cuerpo procede?. El truco para acotar el vector de magnetización transversal a un plano estrecho del cuerpo consiste en no someter al paciente a un campo magnético homogéneo, ya que en ese caso todos los H + precesarían a la misma frecuencia y en consecuencia todos serían alterados por la misma RF. Si, en cambio, el paciente se somete a un gradiente de campo, cada tramo de su cuerpo tendrá los H + precesando a una frecuencia característica, dada por la ecuación de Larmor. Así, variando la RF obtenemos señales de distintos planos del cuerpo, tanto más finos cuanta mayor sea la graduación del 0 B r . El siguiente paso es cortar la radiofrecuencia, recordamos que mandamos pulsos de radiofrecuencia (PRF). Ahora los H + no tienen ninguna energía que les haga girar en fase y poco a poco vuelven a desfasarse, además de recuperar su estado de menor energía. Comienza a crecer la Conceptos básicos de Resonancia Magnética Nuclear 19 magnetización longitudinal a medida que decrece la transversal. La velocidad con la que se producen estos fenómenos (T1 y T2) es característica para tejido. El ordenador es capaz de dar un código de grises para cada intensidad de señal, marcando así cada uno de los tejidos. 2.4 El concepto de vóxel 2.4 El concepto de vóxel2.4 El concepto de vóxel 2.4 El concepto de vóxel El vóxel (Volumetric Pixel o más correctamente Volumetric Picture Element) es una unidad cúbica que compone un objeto tridimensional. Constituye la unidad mínima procesable de una matriz tridimensional y es, por tanto, el equivalente al píxel en un objeto 2D. Figura 2.6. Figura 2.6.Figura 2.6. Figura 2.6. Representación de vóxeles en una matriz tridimensional. Para crear una imagen en tres dimensiones, los vóxeles tienen que sufrir una transformación de opacidad. Esta información da diferentes valores de opacidad a cada vóxel. Esto es importante cuando se han de mostrar detalles interiores de una imagen que quedaría tapada por la capa exterior más opaca de los vóxeles. Las imágenes con vóxeles se usan generalmente en el campo de la medicina y se aplican, por ejemplo, en la tomografía axial computarizada o para las resonancias magnéticas. De este modo, los profesionales pueden obtener un modelo preciso en tres dimensiones del cuerpo humano. Conceptos básicos de Resonancia Magnética Nuclear 20 Figura 2.7. Figura 2.7.Figura 2.7. Figura 2.7. Vóxel en el que se pueden diferenciar distintos tejidos. Actualmente, su uso ya se ha extendido en multitud de campos como la medicina, ingeniería, cine, videojuegos... Al igual que los píxeles, los vóxeles no contienen su posición (x,y,z) en el espacio 3D, sino que esta se deduce por la posición del vóxel dentro del archivo de datos, y tienen un tamaño concreto. 3. 3.3. 3. Material y métodos Material y métodosMaterial y métodos Material y métodos utilizados utilizados utilizados utilizados 3.1 Imágenes médicas 3.1 Imágenes médicas3.1 Imágenes médicas 3.1 Imágenes médicas El desarrollo de la aplicación está enfocada al análisis y manipulación de imágenes médicas correspondientes al encéfalo. Para la prueba del software se han utilizado imágenes anónimas de pacientes de varias edades, con orientación mayormente axial. A la hora de tomar las imágenes, no existe un convenio sobre qué orientación deben tener, por ello el software desarrollado permite reorientar las imágenes a una orientación axial una vez cargadas, previa selección del usuario sobre qué orientación tenía inicialmente. Material y métodos utilizados 21 Figura 3.1. Figura 3.1.Figura 3.1. Figura 3.1. Corte axial (izquierda). Corte coronal (centro). Corte sagital (derecha). Es necesario realizar una reorientación en el caso de que no sean axiales para poder aplicar los métodos desarrollados por el grupo de imagen médica. Las imágenes utilizadas tendrán formato .img y .nii . 3.2 Matlab 3.2 Matlab3.2 Matlab 3.2 Matlab MATLAB (MATrix LABoratory) es un software matemático que ofrece un entorno de desarrollo integrado (IDE) con un lenguaje de programación propio (scripts y ficheros .mex escritos en lenguaje C). Figura 3.2. Figura 3.2.Figura 3.2. Figura 3.2. Pantalla principal de MATLAB. Material y métodos utilizados 22 Entre sus prestaciones básicas se hallan: • Manipulación de matrices. • Representación de datos y funciones. • Implementación de algoritmos • Creación de interfaces de usuario • Comunicación con programas en otros lenguajes y con otros dispositivos hardware. El paquete MATLAB dispone de una herramienta adicional que expande sus prestaciones, y que nos va a resultar de gran ayuda a la hora de crear interfaces de usuario:GUIDE. Figura 3.3. Figura 3.3.Figura 3.3. Figura 3.3. Pantalla de la herramienta GUIDE. El principal motivo por el cual se ha elegido MATLAB como herramienta es precisamente su capacidad y eficacia para el tratamiento de matrices. Las imágenes médicas están almacenadas internamente en forma matricial, y es por lo tanto una herramienta óptima. Otra gran cualidad es la posibilidad de crear scripts, y son precisamente estos los que hemos utilizado a la hora de desarrollar el software. Así mismo, los módulos desarrollados por el grupo de imagen médica están escritos en matlab (los Material y métodos utilizados 23 scripts de matlab tienen extensión .m) junto con, en algunos casos, en lenguaje C (ficheros .mex). 3.3 3.33.3 3.3 SPM SPM SPM SPM Los recientes avances en técnicas de neuroimagen han contribuido a la proliferación de investigaciones y estudios relacionados con este campo. Los más vanguardistas utilizan una herramienta conocida como Statistical Parametric Mapping (SPM). SPM se ha convertido en el estándar de facto para la aplicación de métodos estadísticos en el análisis de imágenes. La herramienta se compone de una colección de funciones y una interfaz gráfica que se ejecuta directamente sobre MATLAB. Figura 3.4. Figura 3.4.Figura 3.4. Figura 3.4. Interfaz gráfica de SPM8. En concreto, la versión que se ha utilizado ha sido SPM8. SPM8 soporta el tratamiento de imágenes de RM, PET y SPECT. Una de las principales ventajas que tiene este software deriva de su carácter open source, ya que son muchos los investigadores (entre los que se incluyen informáticos, médicos, matemáticos, estadísticos, farmacéuticos y físicos) que a nivel mundial trabajan con este sistema. La aportación de cada uno de estos usuarios conlleva un desarrollo constante de la aplicación e incrementa su robustez frente a posibles errores o bugs. Además, SPM8 posee una buena granularidad de funciones que hace posible la creación de Material y métodos utilizados 24 procesos automatizados encapsulados en procesos por lotes (también conocidos como procesos batch o scripts). Es precisamente esta herramienta la que vamos a utilizar para poder leer los ficheros de imagen y obtener la matriz correspondiente para poder aplicar los tratamientos oportunos. Aunque los módulos desarrollados por el grupo de imagen médica hacen uso de más llamadas, para cargar los casos y sacar la información de la matriz utilizaremos principalmente vamos a utilizar estas dos llamadas: spm_vol spm_volspm_vol spm_vol V = spm_vol(P) P - a matrix of filenames. V - a vector of structures containing image volume information. The elements of the structures are: V.fname - the filename of the image. V.dim - the x, y and z dimensions of the volume V.dt - A 1x2 array. First element is datatype (see spm_type). The second is 1 or 0 depending on the endian-ness. V.mat - a 4x4 affine transformation matrix mapping from voxel coordinates to real world coordinates. V.pinfo - plane info for each plane of the volume. V.pinfo(1,:) - scale for each plane V.pinfo(2,:) - offset for each plane The true voxel intensities of the jth image are given by: val*V.pinfo(1,j) + V.pinfo(2,j) V.pinfo(3,:) - offset into image (in bytes). If the size of pinfo is 3x1, then the volume is assumed to be contiguous and each plane has the same scalefactor and offset. spm_read_vols spm_read_volsspm_read_vols spm_read_vols Y = spm_read_vols(V) V - vector of mapped image volumes to read in (from spm_vol) Y - 4D matrix of image data, fourth dimension indexes images Básicamente estas dos llamadas nos van a servir para cargar una imagen a partir de su nombre de fichero ( spm_vol ) y posteriormente utilizar esta variable para cargar los datos propios del volumen ( spm_read_vols ). Con la llamada spm_vol además de cargar la imagen, nos devolverá información como el nombre de la imagen, las dimensiones del volumen y el tamaño de los vóxeles. Material y métodos utilizados 25 3.4 Módulos 3.4 Módulos3.4 Módulos 3.4 Módulos El objetivo final de este proyecto es el desarrollo de una interfaz gráfica de usuario en entorno MATLAB para la integración de los módulos desarrollados por el grupo de imagen médica. En este punto se hará una breve descripción de cada uno de los módulos, dejando una descripción más detallada en los artículos citados en la bibliografía. Filtrado de ruido En el proceso de adquisición de imágenes es prácticamente imposible que se obtengan sin ningún tipo de ruido. El ruido en las adquisiciones puede deberse a varios factores que escapan al factor humano. Algunas de estas causas pueden ser la inhomogeneidades del campo magnético, radiaciones de fondo y movimientos del propio paciente dentro de la máquina. (a) (b) (c) Figura Figura Figura Figura 3.5. 3.5.3.5. 3.5. Imagen sin ruido(a), imagen con un ruido estimado del 9%(b) y los residuos correspondientes (c). Es necesario aplicar un buen filtrado a las imágenes, ya que la clasificación del tejido se hace en base al brillo de este, y el ruido puede causar una confusión a la hora de clasificar y por ello producir resultados indeseados. El módulo de filtrado puede realizar el filtrado propiamente dicho o solamente una estimación del ruido, realizándose esta como paso previo para realizar el filtrado, guardando además del fichero con el volumen sin ruido, otro fichero con los residuos. Resultados 32 En la parte inferior derecha está la zona de herramientas, en la que se puede aplicar el zoom, mover la imagen y, una vez activado el editor, corregir la máscara con una etiqueta y un tamaño de dibujado concreto. Para cambiar el corte, con hacer click izquierdo en una vista del volumen se actualizarán las vistas de los otros dos cortes, de la misma forma, con la rueda del ratón se pueden cambiar los cortes de la vista actual. La otra opción, Batch Process , muestra una interfaz muy simple e intuitiva: Figura 5.3. Figura 5.3.Figura 5.3. Figura 5.3. Pantalla de la opción Batch Process del programa. En la pantalla de esta opción hay tres bloques; para especificar la orientación que tienen todos los volúmenes, seleccionar qué procesos queremos ejecutar y las opciones de guardado. En estas últimas se dan dos posibles casos: • Guardar dentro del subdirectorio results • Seleccionar un directorio por el usuario Resultados 33 En el primer caso se guardarán los resultados (ficheros intermedios e información de la segmentación) dentro del subdirectorio results bajo otro subdirectorio que tendrá como nombre el que se introduzca en el campo de texto. En el segundo caso los resultados se guardarán en la ubicación que especifique el usuario mediante un cuadro de diálogo. 6. 6. 6. 6. Discusión y conclusiones Discusión y conclusionesDiscusión y conclusiones Discusión y conclusiones El proyecto comenzó a mediados de noviembre con la lectura de documentación sobre matlab y la programación de scripts para este, para posteriormente empezar a desarrollar una primera aproximación a la interfaz, que ha llevado unos ocho meses de desarrollo en etapas más o menos activas. Como resultado de este proyecto se ha obtenido un software usable, en el que se ha hecho todo lo posible por automatizar todas las tareas de tratamiento de imágenes médicas, y hacerlo lo más sencillo posible. Durante el desarrollo ha habido ciertas partes que han resultado más problemáticas, como la redimensión de los volúmenes independientemente de su tamaño original o las opciones de pincel y la iniciación a matlab. Hasta ahora el uso de los métodos desarrollados por el grupo de imagen médica sólo se podían utilizar mediante la línea de comandos de matlab, por lo que se espera que esta herramienta resulte de gran utilidad para los operarios. El software está abierto a futuras mejoras y actualizaciones en base a las propias mejoras en los métodos que integra. Personalmente, creo que ha sido un proyecto bonito e interesante, en el que he aprendido muchas cosas sobre imagen médica y matlab, y que espero sirva de gran utilidad. Agradecimientos 34 7. Agradecimientos 7. Agradecimientos7. Agradecimientos 7. Agradecimientos Esta es sin duda la parte más personal del documento, y en la que más nos pensamos qué decir y cómo decirlo, y a la vez quizás la parte más complicada. Sin duda alguna, a las personas que más tengo que agradecer en la vida son a mi madre y a mi abuela. Siempre han estado ahí, desde pequeño, apoyándome y dándome consejos, educándome como mejor han creído y sacándome adelante siempre pese a las circunstancias. Si no fuese por ellas jamás podría haber estudiado lo que me gusta, ni estar donde estoy ahora, y ha sido por ellas por las que muchas veces he sacado fuerzas cuando pensaba que algo era demasiado complicado, algo que necesitaba esfuerzo y sacrificio. Han sabido enseñarme cómo deben ser las cosas y cómo no deben serlo, a luchar por un futuro y ser alguien. Les debo la vida. A parte de mi familia, mi padre y mi tío, que también me han aconsejado, me han apoyado y han confiado en mí, y cada uno a su manera me han enseñado hacia dónde debía dirigir mi futuro. A mi novia, Sylvia, que me ha apoyado siempre, ha sabido escucharme, me ha aconsejado y ha estado a mi lado en todo momento, siempre que lo he necesitado, y que es una parte fundamental en mi vida. Te quiero. Y sin duda a mi director de proyecto, José, que ha hecho lo imposible por ayudarme siempre que me he atascado con algo, mostrando interés por que progresase. Gracias a él he aprendido muchas cosas con este proyecto, tanto en el ámbito informático como en el médico, y me ha dedicado un montón de horas de trabajo. Si todos los directores fuesen así, tendríamos un bonito recuerdo de nuestro último escalón al finalizar la carrera en lugar de un mero trámite para terminarla. Gracias José. Y a mi abuelo, al que siempre llevaré en mi memoria. Bibliografía y referencias 35 8. Bibliografía y referencias 8. Bibliografía y referencias8. Bibliografía y referencias 8. Bibliografía y referencias Artículos científicos José V.Manjón, Juan J. Lull, José Carbonell-Caballero, Gracían GarcíaMartí, Luis Martí-Bonmatí, Montserrat Robles. 2007. A nonparametricMRI inhomogeneity correction method. Medical image analysis 11: 336-345. Pierrick Coupé, Pierre Yger, Sylvain Prima, Pierre Hellier, Charles Kervrann, Christian Barillot. 2008. An optimized blockwise nonlocal means denoising filter for 3-D magnetic resonante images. José V.Manjón, Juan J. Lull, José Carbonell-Caballero, Gracían GarcíaMartí, Luis Martí-Bonmatí, Montserrat Robles. 2008. MRI denoising using non-local means. Medical image analysis 12: 514-523. José V.Manjón, Jussi Tohka, Gracián García-Martí, José CarbonellCaballero, Juan J.Lull, Luís Matí-Bonmatí, Montserrat Robles. 2008. Robust MRI brain tissue parameter estimation by multistage outlier rejection. Magnetic resonance in medicine 59: 866-873. Tutoriales Mª Victoria Lapuerta González, Ana Laverón Simavilla. Introducción a MATLAB. Gonzalo Fernández de Córdoba Martos. Creación de interfaces gráficas de usuario con MATLAB. Ashok Krishnamurthy, Siddharth Samsi. MATLAB for engineering applications. Diego Orlando Barragán Guerrero. Manual de interfaz gráfica de usuario en MATLAB. Emma Díez. Imágenes de resonancia magnética nuclear. ANEXO 36 ANEXO ANEXOANEXO ANEXO En esta sección se mostrará el código de los siguientes scripts implementados para el desarrollo del software: • brainseg.m • gui.m • batch.m • loadFile.m • reorientation.m • volReformat.m • resizeVol.m • save_volume.m • run_batch.m Brainseg.m Brainseg.mBrainseg.m Brainseg.m function varargout = brainseg(varargin) % Begin initialization code – DO NOT EDIT if(exist('defaults.mat','file')) load defaults; addpath(spmpath); addpath(vbmpath); end; isspm = findstr(path,'spm8'); if(isempty(isspm)) spmpath=uigetdir(pwd,'Select SPM8 directory'); addpath(spmpath); end isvbm = findstr(path,'vbm8'); if(isempty(isvbm)) vbmpath=uigetdir(pwd,'Select VBM8 directory'); addpath(vbmpath); save defaults spmpath vbmpath end % Add additional subdirectories actual=pwd; addpath([actual, '\modules']); addpath([actual, '\modules\Denoising']); addpath([actual, '\modules\BrainExtraction']); addpath([actual, '\modules\Clasification']); % Begin initialization code - DO NOT EDIT gui_Singleton = 1; gui_State = struct('gui_Name', mfilename, ... 'gui_Singleton', gui_Singleton, ... ANEXO 37 'gui_OpeningFcn', @brainseg_OpeningFcn, ... 'gui_OutputFcn', @brainseg_OutputFcn, ... 'gui_LayoutFcn', [] , ... 'gui_Callback', []); if nargin && ischar(varargin{1}) gui_State.gui_Callback = str2func(varargin{1}); end if nargout [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:}); else gui_mainfcn(gui_State, varargin{:}); end % End initialization code - DO NOT EDIT % --- Executes just before brainseg is made visible. function brainseg_OpeningFcn(hObject, eventdata, handles, varargin) % Choose default command line output for brainseg handles.output = hObject; % Update handles structure guidata(hObject, handles); % UIWAIT makes brainseg wait for user response (see UIRESUME) % uiwait(handles.figure1); global arg arg.img_name = []; arg.axial=0; arg.coronal=0; arg.sagittal=0; arg.s=[]; arg.csag=1; arg.ccor=1; arg.cax=1; axes(handles.axesLogo) handles.imagen=imread('logo.tif'); imagesc(handles.imagen) axis off %End of global var args % --- Outputs from this function are returned to the command line. function varargout = brainseg_OutputFcn(hObject, eventdata, handles) % Get default command line output from handles structure varargout{1} = handles.output; % --- Executes on button press in guiButton. function guiButton_Callback(hObject, eventdata, handles) close; gui; ANEXO 38 % --- Executes on button press in batchButton. function batchButton_Callback(hObject, eventdata, handles) close; batch; gui.m gui.mgui.m gui.m function varargout = gui(varargin) % Begin initialization code - DO NOT EDIT gui_Singleton = 1; gui_State = struct('gui_Name', mfilename, ... 'gui_Singleton', gui_Singleton, ... 'gui_OpeningFcn', @gui_OpeningFcn, ... 'gui_OutputFcn', @gui_OutputFcn, ... 'gui_LayoutFcn', [] , ... 'gui_Callback', []); if nargin && ischar(varargin{1}) gui_State.gui_Callback = str2func(varargin{1}); end if nargout [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:}); else gui_mainfcn(gui_State, varargin{:}); end % End initialization code - DO NOT EDIT % --- Executes just before gui is made visible. function gui_OpeningFcn(hObject, eventdata, handles, varargin) % Choose default command line output for gui handles.output = hObject; % Update handles structure guidata(hObject, handles); % UIWAIT makes gui wait for user response (see UIRESUME) % uiwait(handles.guiWindow); global vol; global V; global arg; global res; global f; global origin; global volImg; global volOriginal; global cutAxial; global cutSagittal; global cutCoronal; global left; global rect; ANEXO 39 global brillo; global sizeVolImg; global maxSize; global overlay; global map; global enable; global resultspath; global fname; global limit_a; global limit_s; global limit_c; global mapa; global right; global pencil; % Initialize tools section set(handles.disable_editor_menu, 'Enable', 'off'); % Both menus are disabled by default set(handles.save_mask_menu, 'Enable', 'off'); % End of init. enable=0; right=0; left=0; brillo=0; overlay=0; pencil=0; if(isempty(arg.img_name)) volOriginal=zeros(64,64,64); res=[1,1,1]; arg.axial = 1; %Disable tool menus if acquisition is not loaded set(handles.denoising_menu, 'Enable', 'off'); set(handles.correction_menu, 'Enable', 'off'); set(handles.extraction_menu, 'Enable', 'off'); set(handles.classification_menu, 'Enable', 'off'); set(handles.segment_menu, 'Enable', 'off'); set(handles.save_options_menu, 'Enable', 'off'); set(handles.buttonZoom, 'Enable', 'off'); set(handles.buttonPan, 'Enable', 'off'); set(handles.buttonPencil, 'Enable', 'off'); set(handles.popupmenu_size, 'Enable', 'off'); set(handles.popupmenu_label, 'Enable', 'off'); %Disable overlay options set(handles.load_overlay_menu, 'Enable', 'off'); set(handles.unload_overlay_menu, 'Enable', 'off'); set(handles.transparency_menu, 'Enable', 'off'); else V = spm_vol(arg.img_name); volOriginal=spm_read_vols(V); res(1)=abs(V.mat(1,1)); res(2)=abs(V.mat(2,2)); ANEXO 40 res(3)=abs(V.mat(3,3)); origin(1)=V.mat(1,4); origin(2)=V.mat(2,4); origin(3)=V.mat(3,4); set(handles.denoising_menu, 'Enable', 'on'); set(handles.correction_menu, 'Enable', 'on'); set(handles.extraction_menu, 'Enable', 'on'); set(handles.classification_menu, 'Enable', 'on'); set(handles.segment_menu, 'Enable', 'on'); set(handles.save_options_menu, 'Enable', 'on'); set(handles.buttonZoom, 'Enable', 'on'); set(handles.buttonPan, 'Enable', 'on'); %Enable overlay options set(handles.load_overlay_menu, 'Enable', 'on'); set(handles.unload_overlay_menu, 'Enable', 'on'); set(handles.transparency_menu, 'Enable', 'on'); %Set saving as axial by default set(handles.save_axial_menu, 'Checked', 'on'); % Set information data in information box [infopath,infoname,infoext]=fileparts(arg.img_name); infoFilename=[infoname,infoext]; set(handles.text_info, 'String', { ['Filename: ', infoFilename ] }); set(handles.text_info, 'String', [get(handles.text_info, 'String'), {['Volume Size: ', num2str(size(volOriginal)) ]} ]); end stemp=size(volOriginal); % If the acquisition is axial if (arg.axial == 1) colormap(gray); arg.s=size(volOriginal); % If all voxels are the same size, then call resizeVol function if (res(1) == res(2) && res(3)==res(1)) volImg=resizeVol(volOriginal); else resMin=min(res); % Store the minimum value of res f=res/resMin; % Calculate new resolution value and call to the interpolation function volTemp=volReformat(volOriginal, f); volImg=resizeVol(volTemp); end % createcolormap Imin=round(min(volImg(:))); ANEXO 41 volImg=volImg-Imin; Imax=round(max(volImg(:))); volImg=volImg*255/Imax; Imax=255; mapa=zeros(Imax,3); % If adquisition is axial copy Original volume in new vol vol=volOriginal; for i=1:Imax mapa(i,1)=(i-1)/(Imax-1); mapa(i,2)=(i-1)/(Imax-1); mapa(i,3)=(i-1)/(Imax-1); end colormap(mapa); sizeVolImg=size(volImg); % For each image we'll show the central cut without labels axes(handles.axesAxial); cutAxial=round((sizeVolImg(3)/2)); image(imrotate(volImg(:,:,cutAxial),90)); set(handles.axesAxial, 'Tag', 'axial'); axis off; axes(handles.axesSagital); cutSagittal=round((sizeVolImg(1)/2)); sag(:,:)=volImg(cutSagittal,:,:); image(imrotate(sag,90)); set(handles.axesSagital, 'Tag', 'sagittal'); axis off; axes(handles.axesCoronal); cutCoronal=round((sizeVolImg(2)/2)); cor(:,:)=volImg(:,cutCoronal,:); image(imrotate(cor,90)); set(handles.axesCoronal, 'Tag', 'coronal'); axis off; end % If the acquisition is coronal if (arg.coronal == 1) colormap(gray); %Create new matrix to store fixed acquisition in axial orientation vol=zeros(stemp(1),stemp(3),stemp(2)); % Copy axials images in the new matrix for i=1:stemp(2) vol(:,:,i)=imrotate(squeeze(volOriginal(:,i,:)),180); end tt=res(2); res(2)=res(3); res(3)=tt; ANEXO 48 % --- Executes on mouse press over figure background, over a disabled or % --- inactive control, or over an axes background. function guiWindow_WindowButtonDownFcn(hObject, eventdata, handles) global volImg; global volImgFusion; global sizeVolImg; global left; global f; global rect; global overImg; global volImgOverlay; global offsetx; global offsety; global offsetz; global cutAxial; global cutSagittal; global cutCoronal; global arg; global enable; global pencil; global right; global sizeOverImg; axedSelected=get(gca, 'Tag'); pos=get(gca,'CurrentPoint'); x=round(pos(1,1)); y=round(pos(1,2)); if strcmp(get(gcf,'SelectionType'),'normal') left=1; right=0; end if strcmp(get(gcf,'SelectionType'),'alt') left=0; right=1; end % If editor is enable and we want to draw if(enable == 1 && pencil==1) if ( ((x>0) && (x<sizeVolImg(1))) && ((y>0) && (y<sizeVolImg(1))) ) switch axedSelected case 'axial' cutCoronal=(round((sizeVolImg(2)-y+1))); cutSagittal=round(((x))); case 'coronal' cutAxial=(sizeVolImg(2)-y)+1; cutSagittal=round(((x))); case 'sagittal' ANEXO 49 cutAxial=(sizeVolImg(2)-y+1); cutCoronal=(round((x))); end x1=round((cutSagittal-offsetx)/f(1)); y1=round((cutCoronal-offsety)/f(2)); z1=round((cutAxial-offsetz)/f(3)); sizeSelected = get(handles.popupmenu_size, 'Value')-1; % Left click if (left==1) labelMax=max(overImg(:)); labelValue=get(handles.popupmenu_label, 'Value'); for i=-sizeSelected:sizeSelected for j=-sizeSelected:sizeSelected for k=-sizeSelected:sizeSelected newX=x1+i; newY=y1+j; newZ=z1+k; d = sqrt(i*i + j*j + k*k); if ( d > sizeSelected) continue; end if( (newX>0 && newX<=sizeOverImg(1)) && (newY>0 && newY<=sizeOverImg(2)) && (newZ>0 && newZ<=sizeOverImg(3)) ) overImg(newX, newY, newZ)=labelValue; end newSagittal=cutSagittal+i; newCoronal=cutCoronal+j; newAxial=cutAxial+k; if ( ((newSagittal>0) && (newSagittal<sizeVolImg(1))) && ((newCoronal>0) && (newCoronal<sizeVolImg(1))) && ((newAxial>0) && (newAxial<sizeVolImg(1)))) volImgOverlay(newSagittal, newCoronal, newAxial)=(255*labelValue)/labelMax; volImgFusion(newSagittal, newCoronal, newAxial)=129 + (volImg(newSagittal, newCoronal, newAxial).*(volImgOverlay(newSagittal, newCoronal, newAxial)))/(256); end end end end ANEXO 50 if (labelValue > labelMax) change_alfa(labelValue); end end % End of if left % Right click if(right==1) for i=-sizeSelected:sizeSelected for j=-sizeSelected:sizeSelected for k=-sizeSelected:sizeSelected newX=x1+i; newY=y1+j; newZ=z1+k; d = sqrt(i*i + j*j + k*k); if ( d > sizeSelected) continue; end if( (newX>0 && newX<=sizeOverImg(1)) && (newY>0 && newY<=sizeOverImg(2)) && (newZ>0 && newZ<=sizeOverImg(3)) ) overImg(newX, newY, newZ)=0; end newSagittal=cutSagittal+i; newCoronal=cutCoronal+j; newAxial=cutAxial+k; volImgOverlay(newSagittal, newCoronal, newAxial)=0; volImgFusion(newSagittal, newCoronal, newAxial)=0.5*volImg(newSagittal, newCoronal, newAxial); end end end end % End of if right axes(handles.axesAxial); image(imrotate(volImgFusion(:,:,cutAxial),90)); set(handles.axesAxial, 'Tag', 'axial'); axis off; zoom(1); axes(handles.axesSagital); sag(:,:)=volImgFusion(cutSagittal,:,:); image(imrotate(sag,90)); set(handles.axesSagital, 'Tag', 'sagittal'); axis off; zoom(1); axes(handles.axesCoronal); cor(:,:)=volImgFusion(:,cutCoronal,:); image(imrotate(cor,90)); set(handles.axesCoronal, 'Tag', 'coronal'); axis off; ANEXO 51 zoom(1); end end if pencil==0 && strcmp(get(gcf,'SelectionType'),'alt') if(strcmp(axedSelected,'axial')) rect=getrect(handles.axesAxial); rect=round(rect); reg=volImg(rect(1):(rect(1)+rect(3)),rect(2):(rect(2)+rect(4)),cutAxia l); m1=min(reg(:)); m2=max(reg(:)); CambiarBrillo(m1,m2); return; end if(strcmp(axedSelected,'sagittal')) rect=getrect(handles.axesSagital); rect=round(rect); reg=volImg(rect(1):(rect(1)+rect(3)),rect(2):(rect(2)+rect(4)),cutSagi ttal); m1=min(reg(:)); m2=max(reg(:)); CambiarBrillo(m1,m2); return; end if(strcmp(axedSelected,'coronal')) rect=getrect(handles.axesCoronal); rect=round(rect); reg=volImg(rect(1):(rect(1)+rect(3)),rect(2):(rect(2)+rect(4)),cutCoro nal); m1=min(reg(:)); m2=max(reg(:)); CambiarBrillo(m1,m2); return; end end position_update(handles, x, y, axedSelected); if(length(arg.img_name) ~= 0) switch axedSelected case 'axial' axes(handles.axesAxial); set(handles.textCoordZ, 'String', cutAxial-1); case 'coronal' axes(handles.axesCoronal); set(handles.textCoordY, 'String', cutCoronal-1); case 'sagittal' axes(handles.axesSagital); set(handles.textCoordX, 'String', cutSagittal-1); end end ANEXO 52 % -------------------------------------------------------------------- function new_acquisition_Callback(hObject, eventdata, handles) global overlay; overlay=0; loadFile; % --- Executes on mouse motion over figure - except title and menu. function guiWindow_WindowButtonMotionFcn(hObject, eventdata, handles) global volImg; global volImgFusion; global volImgOverlay; global overImg; global left; global right; global cutAxial; global cutSagittal; global cutCoronal; global arg; global f; global offsetx; global offsety; global offsetz; global enable; global pencil; global sizeVolImg; global sizeOverImg; s=sizeVolImg; axedSelected=get(gca, 'Tag'); posPointer=get(gca,'CurrentPoint'); if(length(arg.img_name)~=0) switch axedSelected case 'axial' xP=round(posPointer(1,1)); yP=(s(2)-round(posPointer(1,2))); if(xP>0 && xP<=s(1) && yP>0 && yP<=s(2)) set(handles.textCoordX, 'String', xP-1); set(handles.textCoordY, 'String', yP-1); set(handles.textCoordZ, 'String', cutAxial-1); set(handles.valor, 'String', volImg(xP,yP,cutAxial)); end case 'coronal' xP=round(posPointer(1,1)); yP=(s(3)-round(posPointer(1,2))); if(xP>0 && xP<=s(1) && yP>0 && yP<=s(3)) set(handles.textCoordX, 'String', xP-1); set(handles.textCoordZ, 'String', yP-1); set(handles.textCoordY, 'String', cutCoronal-1); set(handles.valor, 'String', volImg(xP,cutCoronal,yP)); end ANEXO 53 case 'sagittal' xP=round(posPointer(1,1)); yP=(s(3)-round(posPointer(1,2))); if(xP>0 && xP<=s(1) && yP>0 && yP<=s(3)) set(handles.textCoordY, 'String', xP-1); set(handles.textCoordZ, 'String', yP-1); set(handles.textCoordX, 'String', cutSagittal-1); set(handles.valor, 'String', volImg(cutSagittal,xP,yP)); end end end axedSelected=get(gca, 'Tag'); pos=get(gca,'CurrentPoint'); x=round(pos(1,1)); y=round(pos(1,2)); if ( (left == 1 || right==1) && (enable == 1) && (pencil == 1)) if ( ((x>0) && (x<sizeVolImg(1))) && ((y>0) && (y<sizeVolImg(1)))) switch axedSelected case 'axial' cutCoronal=(round((sizeVolImg(2)-y+1))); cutSagittal=round(((x))); case 'coronal' cutAxial=(sizeVolImg(2)-y+1); cutSagittal=round(((x))); case 'sagittal' cutAxial=(sizeVolImg(2)-y+1); cutCoronal=(round((x))); end x1=round((cutSagittal-offsetx)/f(1)); y1=round((cutCoronal-offsety)/f(2)); z1=round((cutAxial-offsetz)/f(3)); sizeSelected = get(handles.popupmenu_size, 'Value')-1; if(left==1) labelMax=max(overImg(:)); labelValue=get(handles.popupmenu_label, 'Value'); for i=-sizeSelected:sizeSelected ANEXO 54 for j=-sizeSelected:sizeSelected for k=-sizeSelected:sizeSelected newX=x1+i; newY=y1+j; newZ=z1+k; d = sqrt(i*i + j*j + k*k); if ( d > sizeSelected) continue; end if( (newX>0 && newX<=sizeOverImg(1)) && (newY>0 && newY<=sizeOverImg(2)) && (newZ>0 && newZ<=sizeOverImg(3)) ) overImg(newX, newY, newZ)=labelValue; end newSagittal=cutSagittal+i; newCoronal=cutCoronal+j; newAxial=cutAxial+k; if ( ((newSagittal>0) && (newSagittal<sizeVolImg(1))) && ((newCoronal>0) && (newCoronal<sizeVolImg(1))) && ((newAxial>0) && (newAxial<sizeVolImg(1)))) volImgOverlay(newSagittal, newCoronal, newAxial)=(255*labelValue)/labelMax; volImgFusion(newSagittal, newCoronal, newAxial)=129 + (volImg(newSagittal, newCoronal, newAxial).*(volImgOverlay(newSagittal, newCoronal, newAxial)))/(256); end end end end if (labelValue > labelMax) change_alfa(labelValue); end end if(right==1) for i=-sizeSelected:sizeSelected for j=-sizeSelected:sizeSelected for k=-sizeSelected:sizeSelected newX=x1+i; newY=y1+j; newZ=z1+k; d = sqrt(i*i + j*j + k*k); if ( d > sizeSelected) continue; end ANEXO 55 if( (newX>0 && newX<=sizeOverImg(1)) && (newY>0 && newY<=sizeOverImg(2)) && (newZ>0 && newZ<=sizeOverImg(3)) ) overImg(newX, newY, newZ)=0; end newSagittal=cutSagittal+i; newCoronal=cutCoronal+j; newAxial=cutAxial+k; volImgOverlay(newSagittal, newCoronal, newAxial)=0; volImgFusion(newSagittal, newCoronal, newAxial)=0.5*volImg(newSagittal, newCoronal, newAxial); end end end end axes(handles.axesAxial); image(imrotate(volImgFusion(:,:,cutAxial),90)); set(handles.axesAxial, 'Tag', 'axial'); axis off; zoom(1); axes(handles.axesSagital); sag(:,:)=volImgFusion(cutSagittal,:,:); image(imrotate(sag,90)); set(handles.axesSagital, 'Tag', 'sagittal'); axis off; zoom(1); axes(handles.axesCoronal); cor(:,:)=volImgFusion(:,cutCoronal,:); image(imrotate(cor,90)); set(handles.axesCoronal, 'Tag', 'coronal'); axis off; zoom(1); end end if(left==1 && pencil == 0) position_update(handles, x, y, axedSelected); switch axedSelected case 'axial' axes(handles.axesAxial); case 'coronal' axes(handles.axesCoronal); case 'sagittal' axes(handles.axesSagital); end end % --- Executes on mouse press over figure background, over a disabled or ANEXO 56 % --- inactive control, or over an axes background. function guiWindow_WindowButtonUpFcn(hObject, eventdata, handles) global left; global right; right=0; left=0; function CambiarBrillo(min,max) global overlay; global mapa; min=ceil(min); if (min == 0) min=1; end max=ceil(max); if (max == 0 || max == 1) max=2; end for i=1:min mapa(i,1)=0; mapa(i,2)=0; mapa(i,3)=0; end for i=min:max mapa(i,1)=(i-min+1)/(max-min+1); mapa(i,2)=(i-min+1)/(max-min+1); mapa(i,3)=(i-min+1)/(max-min+1); end for i=max:256 mapa(i,1)=1; mapa(i,2)=1; mapa(i,3)=1; end if (overlay == 1) change_alfa(0); else colormap(mapa); end function run_denoising_menu_Callback(hObject, eventdata, handles) nV=save_volume; denoise(nV.fname); function noise_estimation_menu_Callback(hObject, eventdata, handles) global volOriginal; ANEXO 57 disp('Estimating noise...'); h=RicianSTD(volOriginal); disp(['Noise estimated: ', num2str(h)]); function correction_menu_Callback(hObject, eventdata, handles) nV=save_volume; NUCorrect(nV.fname); function extraction_menu_Callback(hObject, eventdata, handles) nV=save_volume; brainextract(nV.fname); function classification_menu_Callback(hObject, eventdata, handles) nV=save_volume; clasify(nV.fname); function segment_menu_Callback(hObject, eventdata, handles) nV=save_volume; nV.fname segment(nV.fname); [path, name, ext]=fileparts(nV.fname); fileName=[path, '\crisp_', name, ext]; load_mask(fileName); function showVolume(overlay, handles) global cutAxial; global cutCoronal; global cutSagittal; global volImg; global volImgFusion; if (overlay == 1) axes(handles.axesAxial); image(imrotate(volImgFusion(:,:,cutAxial),90)); set(handles.axesAxial, 'Tag', 'axial'); axis off; axes(handles.axesSagital); sag(:,:)=volImgFusion(cutSagittal,:,:); image(imrotate(sag,90)); set(handles.axesSagital, 'Tag', 'sagittal'); axis off; axes(handles.axesCoronal); cor(:,:)=volImgFusion(:,cutCoronal,:); image(imrotate(cor,90)); set(handles.axesCoronal, 'Tag', 'coronal'); axis off; elseif (overlay == 0) ANEXO 64 end colormap(mcm); showVolume(overlay, handles); else if (enable == 0) enable=1; set(handles.enable_editor_menu, 'Enable', 'off'); set(handles.disable_editor_menu, 'Enable', 'on'); set(handles.save_mask_menu, 'Enable', 'on'); set(handles.buttonPencil, 'Enable', 'on'); set(handles.popupmenu_size, 'Enable', 'on'); set(handles.popupmenu_label, 'Enable', 'on'); end end end function disable_editor_menu_Callback(hObject, eventdata, handles) global enable; if (enable == 1) enable=0; set(handles.enable_editor_menu, 'Enable', 'on'); set(handles.disable_editor_menu, 'Enable', 'off'); set(handles.save_mask_menu, 'Enable', 'off'); set(handles.buttonPencil, 'Enable', 'off'); set(handles.popupmenu_size, 'Enable', 'off'); set(handles.popupmenu_label, 'Enable', 'off'); end function save_axial_menu_Callback(hObject, eventdata, handles) set(handles.save_axial_menu, 'Checked', 'on'); set(handles.save_original_menu, 'Checked', 'off'); function save_original_menu_Callback(hObject, eventdata, handles) set(handles.save_axial_menu, 'Checked', 'off'); set(handles.save_original_menu, 'Checked', 'on'); function exit_menu_Callback(hObject, eventdata, handles) close all function save_mask_menu_Callback(hObject, eventdata, handles) global overImg; global VOverlay; ANEXO 65 [path,filename,ext]=fileparts(VOverlay.fname); filename=[filename,ext]; [fileNameOverlay path]=uiputfile(filename, 'Save your overlay file'); if (fileNameOverlay ~= 0) VOverlay.fname=fileNameOverlay; spm_write_vol(VOverlay, overImg); disp('Mask saved!'); end % --- Executes on button press in buttonZoom. function buttonZoom_Callback(hObject, eventdata, handles) global volImgFusion; global cutCoronal; global cutAxial; global cutSagittal; global limit_a; global limit_s; global limit_c; if ( get(handles.buttonZoom, 'Value') == 1 ) set(handles.buttonPan, 'Value', 0); set(handles.buttonPencil, 'Value', 0); set(handles.buttonPencil, 'ForegroundColor', 'Black'); pan off; limit_a = getappdata(handles.axesAxial,'zoom_zoomOrigAxesLimits'); limit_s = getappdata(handles.axesSagital,'zoom_zoomOrigAxesLimits'); limit_c = getappdata(handles.axesCoronal,'zoom_zoomOrigAxesLimits'); zoom on; zoom reset; else zoom off; end function buttonPan_Callback(hObject, eventdata, handles) if ( get(handles.buttonPan, 'Value') == 1 ) set(handles.buttonZoom, 'Value', 0); set(handles.buttonPencil, 'Value', 0); set(handles.buttonPencil, 'ForegroundColor', 'Black'); zoom off; pan on; else pan off; end function buttonPencil_Callback(hObject, eventdata, handles) ANEXO 66 global pencil; global volImgFusion; global cutCoronal; global cutAxial; global cutSagittal; global limit_a; global limit_s; global limit_c; if ( get(handles.buttonPencil, 'Value') == 1 ) set(handles.buttonPencil, 'ForegroundColor', 'Red'); set(handles.buttonZoom, 'Value', 0); zoom off; set(handles.buttonPan, 'Value', 0); pan off; pencil = 1; zoom reset; else set(handles.buttonPencil, 'ForegroundColor', 'Black'); pencil = 0; setappdata(handles.axesAxial,'zoom_zoomOrigAxesLimits',limit_a); setappdata(handles.axesSagital,'zoom_zoomOrigAxesLimits',limit_s); setappdata(handles.axesCoronal,'zoom_zoomOrigAxesLimits',limit_c); axes(handles.axesAxial); image(imrotate(volImgFusion(:,:,cutAxial),90)); set(handles.axesAxial, 'Tag', 'axial'); axis off; axes(handles.axesSagital); sag(:,:)=volImgFusion(cutSagittal,:,:); image(imrotate(sag,90)); set(handles.axesSagital, 'Tag', 'sagittal'); axis off; axes(handles.axesCoronal); cor(:,:)=volImgFusion(:,cutCoronal,:); image(imrotate(cor,90)); set(handles.axesCoronal, 'Tag', 'coronal'); axis off; end % --- Executes during object creation, after setting all properties. function popupmenu_size_CreateFcn(hObject, eventdata, handles) if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) set(hObject,'BackgroundColor','white'); end ANEXO 67 % --- Executes during object creation, after setting all properties. function popupmenu_label_CreateFcn(hObject, eventdata, handles) if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) set(hObject,'BackgroundColor','white'); end function documentation_menu_Callback(hObject, eventdata, handles) dir=[pwd, '\doc\index.html']; open(dir); function about_menu_Callback(hObject, eventdata, handles) open about.fig; batch.m batch.mbatch.m batch.m function varargout = batch(varargin) % Begin initialization code - DO NOT EDIT gui_Singleton = 1; gui_State = struct('gui_Name', mfilename, ... 'gui_Singleton', gui_Singleton, ... 'gui_OpeningFcn', @batch_OpeningFcn, ... 'gui_OutputFcn', @batch_OutputFcn, ... 'gui_LayoutFcn', [] , ... 'gui_Callback', []); if nargin && ischar(varargin{1}) gui_State.gui_Callback = str2func(varargin{1}); end if nargout [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:}); else gui_mainfcn(gui_State, varargin{:}); end % End initialization code - DO NOT EDIT % --- Executes just before batch is made visible. function batch_OpeningFcn(hObject, eventdata, handles, varargin) % Choose default command line output for batch handles.output = hObject; % Update handles structure guidata(hObject, handles); global checkboxActive; checkboxActive = 0; global numFiles; numFiles=0; global opts; ANEXO 68 opts.segment=0; opts.denoising=0; opts.correction=0; opts.extraction=0; opts.clasify=0; opts.remove=0; opts.reorient=1; opts.orientation='axial'; % --- Outputs from this function are returned to the command line. function varargout = batch_OutputFcn(hObject, eventdata, handles) % Get default command line output from handles structure varargout{1} = handles.output; function select_files_batch_menu_Callback(hObject, eventdata, handles) global numFiles; global files; [fileName pathBatch]=uigetfile({'*.img; *.nii'}, 'Choose your acquisition file(s)', 'MultiSelect', 'on'); tam=size(fileName); if(tam(1,2) >1) if ( ischar(fileName(1,1)) ) numFiles=1; set(handles.text_files_loaded, 'String', [num2str(numFiles), ' file loaded']); files=[pathBatch,char(fileName)]; else selected=size(fileName); numFiles=selected(1,2); set(handles.text_files_loaded, 'String', [num2str(numFiles), ' files loaded']); for i=1:numFiles filestemp(i,:)=([pathBatch,char(fileName(i))]); end files=filestemp; end end function select_directory_batch_menu_Callback(hObject, eventdata, handles) global files; global numFiles; global opts; directoryPath = uigetdir('Choose your acquisitions folder'); % If we've selected one directory.. ANEXO 69 if (directoryPath ~=0) lista1=ls([directoryPath, '/*.nii']); lista2=ls([directoryPath, '/*.img']); s1=size(lista1); s2=size(lista2); files=char(zeros((s1(1)+s2(1)),max(s1(2),s2(2)))); s=size(files); for i=1:s1(1) files(i,1:length(lista1(i,:)))=lista1(i,:); end for i=(s1(1)+1):s(1) files(i,1:length(lista2(i-s1(1),:)))=lista2(i-s1(1),:); end tam=size(files); numFiles=tam(1,1); if (numFiles == 1) set(handles.text_files_loaded, 'String', [num2str(numFiles), ' file loaded']); elseif (numFiles > 1) set(handles.text_files_loaded, 'String', [num2str(numFiles), ' files loaded']); end for i=1:numFiles filestemp(i,:)=[directoryPath,'\',files(i,:)]; end files=filestemp; end function exit_batch_menu_Callback(hObject, eventdata, handles) global files; global numFiles; clearvar files; clearvar numFiles; close; function editPath_CreateFcn(hObject, eventdata, handles) if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) set(hObject,'BackgroundColor','white'); end function button_set_path_Callback(hObject, eventdata, handles) global savePath; ANEXO 70 savePath = uigetdir('Choose your acquisitions folder'); if (savePath ~= 0) set(handles.editPath, 'String', savePath); end function pushbutton_run_Callback(hObject, eventdata, handles) global checkboxActive; global numFiles; global files; global opts; global pathSave; if( (isempty(get(handles.editPath, 'String'))) && (isempty(get(handles.editBrainseg, 'String'))) ) errordlg('You must specify the destination path before run batch process.', 'No Path Specified'); elseif (checkboxActive == 0) errordlg('You must choose at least one process to run.', 'No Processes Chosen'); elseif (numFiles == 0) errordlg('Load any file before run.', 'No File(s) Loaded'); else warning('OFF', 'MATLAB:MKDIR:DirectoryExists'); % Disable warning disp('Creating directories and storing data in axial orientation...'); if (~isempty(get(handles.editPath, 'String'))) % If user entered a complete path pathSave=get(handles.editPath, 'String'); else userSelect=get(handles.editBrainseg, 'String'); %If user typed a personal dir under brainseg\results\ mkdir(['results\', userSelect]); pathSave=[pwd, '\results\', userSelect]; end opts.nfiles=numFiles; opts.save=pathSave; if (get(handles.radio_axial, 'Value') == 1) opts.orientation='axial'; end if (get(handles.radio_coronal, 'Value') == 1) opts.orientation='coronal'; end if (get(handles.radio_sagittal, 'Value') == 1) opts.orientation='sagittal'; end if (get(handles.checkbox_segment, 'Value') == 1) opts.segment=1; ANEXO 71 end if (get(handles.checkbox_remove_intermedied, 'Value') == 1) opts.remove=1; end if (get(handles.checkbox_denoising, 'Value') == 1) opts.denoising=1; end if (get(handles.checkbox_correction, 'Value') == 1) opts.correction=1; end if (get(handles.checkbox_extraction, 'Value') == 1) opts.extraction=1; end if (get(handles.checkbox_classification, 'Value') == 1) opts.clasify=1; end run_batch(files,opts); close; end function checkbox_denoising_Callback(hObject, eventdata, handles) global checkboxActive; if (get(hObject, 'Value') == 1) checkboxActive=1; else checkboxActive=0; end function checkbox_correction_Callback(hObject, eventdata, handles) global checkboxActive; if (get(hObject, 'Value') == 1) checkboxActive=1; else checkboxActive=0; end function checkbox_extraction_Callback(hObject, eventdata, handles) global checkboxActive; if (get(hObject, 'Value') == 1) checkboxActive=1; else checkboxActive=0; end ANEXO 72 function checkbox_segment_Callback(hObject, eventdata, handles) global checkboxActive; if (get(hObject, 'Value') == 1) set(handles.checkbox_denoising, 'Value', 1); set(handles.checkbox_correction, 'Value', 1); set(handles.checkbox_extraction, 'Value', 1); set(handles.checkbox_classification, 'Value', 1); checkboxActive=1; else set(handles.checkbox_denoising, 'Value', 0); set(handles.checkbox_correction, 'Value', 0); set(handles.checkbox_extraction, 'Value', 0); set(handles.checkbox_classification, 'Value', 0); checkboxActive=0; end function checkbox_classification_Callback(hObject, eventdata, handles) global checkboxActive; if (get(hObject, 'Value') == 1) checkboxActive=1; else checkboxActive=0; end function panel_orientation_SelectionChangeFcn(hObject, eventdata, handles) if (hObject == handles.radio_axial) opts.orientation='axial'; end if (hObject == handles.radio_coronal) opts.orientation='coronal'; end if (hObject == handles.radio_sagittal) opts.orientation='sagittal'; end function panel_orientation_CreateFcn(hObject, eventdata, handles) opts.orientation='axial'; function editBrainseg_CreateFcn(hObject, eventdata, handles) if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) set(hObject,'BackgroundColor','white'); end loadFile.m loadFile.mloadFile.m loadFile.m ANEXO 73 function varargout = loadFile(varargin) % Begin initialization code - DO NOT EDIT gui_Singleton = 1; gui_State = struct('gui_Name', mfilename, ... 'gui_Singleton', gui_Singleton, ... 'gui_OpeningFcn', @loadFile_OpeningFcn, ... 'gui_OutputFcn', @loadFile_OutputFcn, ... 'gui_LayoutFcn', [] , ... 'gui_Callback', []); if nargin && ischar(varargin{1}) gui_State.gui_Callback = str2func(varargin{1}); end if nargout [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:}); else gui_mainfcn(gui_State, varargin{:}); end % End initialization code - DO NOT EDIT function loadFile_OpeningFcn(hObject, eventdata, handles, varargin) % Choose default command line output for loadFile handles.output = hObject; % Update handles structure guidata(hObject, handles); function varargout = loadFile_OutputFcn(hObject, eventdata, handles) % Get default command line output from handles structure varargout{1} = handles.output; function editFilePath_CreateFcn(hObject, eventdata, handles) if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor')) set(hObject,'BackgroundColor','white'); end function btnBrowse_Callback(hObject, eventdata, handles) global arg; [fileName path]=uigetfile({'*.img; *.nii'}, 'Choose your acquisition file'); if ~(fileName == 0) arg.img_name=strcat(path, fileName); V=spm_vol(arg.img_name); vol=spm_read_vols(V); % Show cut 100 preview s=size(vol); axes(handles.axesPreview); cor(:,:)=vol(:,:,round(s(3)/2)); colormap(gray); ANEXO 80 % Check if user has checked inhomogenity correction if (opts.correction == 1) disp(' '); disp('Applying Inhomogenity Correction on file ...'); NUCorrect(axFile); % Update variable path [path,name,ext]=fileparts(axFile); fname=['m_', name, ext]; axFile=[path, '\', fname]; end % Check if user has checked brain extraction if (opts.extraction == 1) disp(' '); disp('Applying Brain Extraction on file ...'); brainextract(axFile); % Update variable path [path,name,ext]=fileparts(axFile); fname=['e_', name, ext]; axFile=[path, '\', fname]; end % Check if user has checked tissue clasificattion if (opts.clasify == 1) disp(' '); disp('Applying Brain Extraction on file ...'); clasify(axFile); end end % Loop end disp(['All data stored under ', opts.save, ' directory.']); warning('ON', 'MATLAB:MKDIR:DirectoryExists'); % Enable again