1. No category

Universidad Politécnica de Madrid
INFORMES:
ITK, MEX y OpenMP
Sandra Rodrı́guez Rodrigo
March 24, 2015
1
1 Índice
1 Índice
2
2 Introducción
3
3 ITK
3.1 Filtro binomial . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.2 Resultado . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3
3
5
4 MEX
4.1 Filtro binomial . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2 Resultado . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
6
6
8
5 OPENMP
5.1 Paralelización bwDistC . . . .
5.2 Resultado . . . . . . . . . . .
5.2.1 Primera comprobación
5.2.2 Segunda comprobación
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
9
9
9
9
10
6 ITKMEX
11
6.1 Filtro binomial . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11
6.2 Resultado . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12
2
2 Introducción
El objetivo de esta práctica es el estudio de las librerı́as ITK (Insight Segmentation and Registration Toolkit) ası́ como de los interfaces MEX y OpenMp para el
tratamiento de imágenes.
3 ITK
3.1 Filtro binomial
Los filtros binomiales son estructuras eficientes basadas en los coeficientes binomiales
resultantes de la implementación de Gauss.
Este filtro usa la convolución gracias a la cual es posible simplificar las matemáticas
y ahorrar tiempo. También resulta muy útil para el procesamiento de imágenes en
algunos casos como para la detención de bordes.
Dado el siguiente código de Matlab:
1
2
3
4
f u n c t i o n imgSal = binomialMatlab ( imgEnt )
maskConv = [ 1 2 1 ; 2 4 2 ; 1 2 1 ] . / 1 6 ;
imgSal = i m f i l t e r ( imgEnt , maskConv , ’ r e p l i c a t e ’ ) ;
\ l a b e l {2}
Se nos pide conseguir una funcionalidad similar usando las librerı́as ITK. El primer
paso será escribir el archivo CMakeLists.txt
1
2
# Indicamos e l nombre que queremos que t e n g a n u e s t r o p r o y e c t o
PROJECT( B i n o m i a l F i l t r o )
3
4
5
6
7
8
9
#Buscamos l i b r e r i a s de ITK
FIND PACKAGE(ITK REQUIRED)
INCLUDE( ${ITK USE FILE } )
# Indicamos que e l p r o y e c t o va a c r e a r una l i b r e r i a e i n d i c a m o s l o s
a r c h i v o s que van a dar l u g a r a l a l i b r e r i a
# Deben e s t a r c r e a d o s a n t e s de e j e c u t a r CMake
ADD EXECUTABLE( Binom
10
11
12
13
14
15
16
17
Binomial . cpp
Binomial . h
)
TARGET LINK LIBRARIES( Binom
$ {ITK LIBRARIES}
)
CMAKE MINIMUM REQUIRED(VERSION 2 . 6 )
3
A partir de ello, se determinan las carpetas donde se encuentran los Fuentes (Binomial.h y Binomial.cpp) y dónde queremos que se guarden los Binarios. Del filtro
solo se configurarán la imagen de entrada y el número de repeticiones que deseamos
que se ejecute. Por lo que, se tendrán como parámetros la imagen cameraman.tif
sobre la que se aplicará, el número de iteraciones elegidas (en nuestro caso 10) y la
imagen filtrada de salida.
1
2
3
4
5
6
7
8
t y p e d e f u n s i g n e d c h a r num ;
// Imagenes en 2D
c o n s t u n s i g n e d i n t numDims = 2 ;
t y p e d e f i t k : : Image< u n s i g n e d char , numDims > ImageEnt ; // Imagen de
Entrada
t y p e d e f i t k : : Image< u n s i g n e d char , numDims > ImageSal ; // Imagen de
Salida
t y p e d e f i t k : : ImageFileReader < ImageEnt > ReaderType ;
t y p e d e f i t k : : I m a g e F i l e W r i t e r < ImageSal > WriterType ;
t y p e d e f i t k : : B i n o m i a l B l u r I m a g e F i l t e r <ImageEnt , ImageSal> B i n B l u r F i l t e r ;
4
3.2 Resultado
Se muestra el resultado con las ITK y se compara con el obtenido tras la implementación desde Matlab con la función binomialMatlab llamada 10 veces.
Figure 3.1: Imagen de entrada
Figure 3.2: ITK - Imagen de salida
Figure 3.3: binomialMatlab - Imagen de salida
El resultado es muy similar pero por diferencias en la implementación, si se observan detenidamente pixel a pixel se verá que no son iguales.
5
4 MEX
Las funciones MEX son muy útiles, ya que nos permiten emplear código escrito en
C/C+ desde Matlab.
4.1 Filtro binomial
Al igual que en el caso anterior, se pide implementar un filtro binomial que tenga la
funcionalidad descrita.
Sin embargo, para llevarlo se tienen en cuenta previamente las siguientes consideraciones:
• Los bordes que conforman el borde será necesario replicar los valores antes de
aplicar el filtro
• Será necesario convertir los valores de unsigned char a coma flotante para el
caso de la entrada, y hacer la conversión opuesta para la salida.
• Debido a que la convolución de la imagen coincide con la correlación de la
misma con la máscara, por ser esta simétrica, se usará correlación a diferencia
de las ITK (además de ser más facil de implementar)
• Se dividirá la imagen en 9 partes para facilitar nuevamente la implementación.
En cuanto al CMakeLists.txt empleado:
1
2
PROJECT( f i l t r o B i n )
IF (WIN32) # Windows
3
4
5
SET(MATLAB LIB [HKEY LOCAL MACHINE\\SOFTWARE\\MathWorks\\MATLAB
\ \ 7 . 1 0 ;MATLABROOT] / e x t e r n / l i b / win64 / m i c r o s o f t )
SET(MATLAB INC ” [HKEY LOCAL MACHINE\\SOFTWARE\\MathWorks\\MATLAB
\ \ 7 . 1 0 ;MATLABROOT] / e x t e r n / i n c l u d e ” )
6
7
8
9
10
11
12
13
14
15
16
17
18
FIND LIBRARY(MATLAB MEX LIBRARY
libmex
$ {MATLAB LIB}
)
FIND LIBRARY(MATLAB MX LIBRARY
libmx
$ {MATLAB LIB}
)
FIND LIBRARY(MATLAB ENG LIBRARY
libeng
$ {MATLAB LIB}
)
19
20
21
22
23
24
25
FIND PATH(MATLAB INCLUDE
”mex . h”
$ {MATLAB INC}
)
INCLUDE DIRECTORIES( ${MATLAB INCLUDE} )
ADD LIBRARY( f i l t B i n MODULE
26
6
27
28
29
30
31
32
u m b r a l i z a c i o n C . cpp
umbralizacionC . h
)
IF (WIN32) # Windows
IF (CMAKE SIZEOF VOID P EQUAL 4 )
# 32 b i t s
SET TARGET PROPERTIES( f i l t B i n PROPERTIES SUFFIX ” . mexw32 ” )
33
34
35
36
37
38
ELSE(CMAKE SIZEOF VOID P EQUAL 4 ) # 64 b i t s
SET TARGET PROPERTIES( f i l t B i n PROPERTIES SUFFIX ” . mexw64 ” )
ENDIF(CMAKE SIZEOF VOID P EQUAL 4 )
ENDIF(WIN32)
39
40
41
42
43
44
45
46
47
# Indicamos que l a l i b r e r i a c r e a d a va a s e r de t i p o MEX
ADD DEFINITIONS(−DMATLAB MEX FILE)
IF (WIN32)
SET TARGET PROPERTIES( f i l t B i n
PROPERTIES
LINK FLAGS ”/ e x p o r t : mexFunction ”
)
ENDIF(WIN32)
48
49
50
51
52
53
54
55
56
57
#I n c o r p o r a r l a s l i b r e r i a s i n t e r f a z e n t r e Matlab y C++
SET(MATLAB LIBRARIES
$ {MATLAB MEX LIBRARY}
$ {MATLAB MX LIBRARY}
$ {MATLAB ENG LIBRARY}
)
TARGET LINK LIBRARIES( f i l t B i n
$ {MATLAB LIBRARIES}
)
Posteriormente, para la obtención del resultado final es necesario implementar la
siguiente función de test:
1
2
3
4
f u n c t i o n imgSal = t e s t ( )
imgEnt=imread ( ’ cameraman . t i f ’ ) ;
imgSal=f i l t B i n ( imgEnt ) ;
imshow ( [ imgEnt , imgSal ] ) ;
5
6
end
Siendo fitlBin el archivo obtenido después de compilar los binarios en Release con
Visual Studio.
En este caso, únicamente se tendrá como parámetro de entrada la imagen a aplicar
el filtro, en nuestro caso la ya tratada cameraman.tif
7
4.2 Resultado
Se obtienen los siguientes resultados:
Figure 4.1: Resultados con MEX
Figure 4.2: Resultados con binomialMatlab
8
5 OPENMP
La interfaz de programación de aplicaciones (API) OpenMp (Open Multi-Processing)
está destinada a la programación multiproceso de memoria compartida en múltiples
plataformas.
Permite el procesamiento en paralelo de programas escritos en C, C++ o Fortran.
Consiste además en un conjunto de directivas para el compilador, unas bibliotecas
de funciones y una serie de variables de entorno que controlan el comportamiento
de la ejecución.
Se basa en la existencia de múltiples hilos de ejecución que comparten una región
de memoria.
5.1 Paralelización bwDistC
Para esta parte, se partirá del código fuente en C++, bwDistC el cual implementa
la función bwdist que calcula para cada pı́xel de entrada la distancia más cercana a
un valor distinto de cero próximo a él. Tanto para la entrada como para la salida,
se usarán matrices lo que facilitará que se ejecute en paralelo.
El objetivo es probar a paralelizar distintas regiones del código y determinar cuál
es la más adecuada mediante una comparación de tiempos.
5.2 Resultado
5.2.1 Primera comprobación
Se hace paralelismo en la inicialización. Para poder llevarlo a cabo será necesario
modificar la siguiente lı́nea de código:
1
2
3
4
5
#pragma omp p a r a l l e l f o r s c h e d u l e ( g u i d e d )
f o r ( i =0; i <numElem ; i ++) {
imgAntes [ i ]= imgIn [ i ] ;
imgOut [ i ] = 0 ;
}
Con la directiva pragma omp for se consigue que las iteraciones del del bucle
for que acompaña sean ejecutadas en paralelo. En el momento de la ejecución, las
iteraciones se distribuyen entre todos los hilos que se lancen.
9
Tras obtener los Binarios, obtener bwDistC.mexw64 compilando desde Release en
Visual Studio y gracias al script scriptBWDist.m con la imagen dada BWHig se ha
obtenido lo siguiente:
El código de C ha tardado una media de: 1.331467 s
El código de C paralelizado ha tardado una media de: 1.272789 s
Se observa que gracias a la paralelización se minoriza el tiempo.
5.2.2 Segunda comprobación
Se realiza nuevamente otra prueba haciendo paralelismo en el bucle for de dilatación:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
#pragma omp p a r a l l e l f o r s c h e d u l e ( g u i d e d )
f o r ( k=0;k<n s l i c e s ; k++) {
f o r ( j =0; j <n c o l ; j ++) {
f o r ( i =0; i <nrow ; i ++) {
i n d = i+j ∗nrow+k∗ numElemSlice ;
i f ( imgDespues [ i n d ]= imgAntes [ i n d ] )
continue ;
i f ( i −1>=0)
i f ( imgDespues [ i n d ]= imgAntes [ ind −1])
continue ;
i f ( i +1<nrow )
i f ( imgDespues [ i n d ]= imgAntes [ i n d +1])
continue ;
i f ( j −1>=0)
i f ( imgDespues [ i n d ]= imgAntes [ ind−nrow ] )
continue ;
i f ( j+1<n c o l )
i f ( imgDespues [ i n d ]= imgAntes [ i n d+nrow ] )
continue ;
i f ( k−1>=0)
i f ( imgDespues [ i n d ]= imgAntes [ ind−numElemSlice ] )
continue ;
i f ( k+1< n s l i c e s )
imgDespues [ i n d ]= imgAntes [ i n d+numElemSlice ] ;
}
}
27
Para este caso se obtiene que el tiempo también es menor:
El código de C paralelizado ha tardado una media de: 1.212789 s
10
6 ITKMEX
El propósito será ver cómo se puede incluir código de las librerı́as ITK en archivos
MEX que serán ejecutados desde MATLAB.
Para ello, a partir del código desarrollado en el primer ejercicio del informe (ITK)
se convertirá en un archivo MEX que posteriormente se probará con el script scriptBinomial.m Además, para este caso, se introducirán manualmente los valores con
los que se quiere trabajar para la umbralización.
6.1 Filtro binomial
El archivo CMakeLists.txt será similar mientras que el archivo de cabecera umbralizacionITKMEX.h contendrá cosas tanto de ITK como de MEX:
1
2
3
4
5
#i n c l u d e
#i n c l u d e
#i n c l u d e
#i n c l u d e
#i n c l u d e
”mex . h”
<i t k I m a g e . h>
<i t k I m a g e R e g i o n I t e r a t o r . h>
<i t k I m a g e R e g i o n C o n s t I t e r a t o r . h>
<i t k B i n a r y T h r e s h o l d I m a g e F i l t e r . h>
Propio de MEX:
1
2
typedef unsigned char PixelEnt ;
typedef unsigned char P i x e l S a l ;
3
Propio de ITK:
1
2
3
c o n s t u n s i g n e d i n t numDims = 2 ;
t y p e d e f i t k : : Image< PixelEnt , numDims > ImageEnt ;
t y p e d e f i t k : : Image< P i x e l S a l , numDims > ImageSal ;
Nuevas clases:
1
2
t y p e d e f i t k : : I m a g e R e g i o n I t e r a t o r < ImageEnt > I t e r a t o r ;
t y p e d e f i t k : : I m a g e R e g i o n C o n s t I t e r a t o r < ImageSal > C o n s t I t e r a t o r ;
Propio de MEX:
1
t y p e d e f i t k : : B i n a r y T h r e s h o l d I m a g e F i l t e r < ImageEnt , ImageSal >
ThresholdFilter ;
2
3
v o i d u m b r a l i z a c i o n ( P i x e l E n t ∗ imgEnt , mwSize numDims , c o n s t mwSize∗
dimImgEnt , P i x e l E n t umbInf , P i x e l E n t umbSup , P i x e l S a l ∗ imgSal ) ;
4
En esta última función sin embargo, se observa un cambio, ya que recibe dos
parámentros más, el valor inferior y superior de umbralización de la imagen de
entrada.
11
6.2 Resultado
Tras seguir los pasos y obtener el archivo umbralizacionITKMEX.mexw64 se desarrolla la siguiente función de test para poder obtener las imágenes:
1
2
3
4
imEnt=imread ( ’ cameraman . t i f ’ ) ;
umbInf =20; umbSup=100;
imSal=umbralizacionITKMEX ( imEnt , umbInf , umbSup ) ;
imshow ( [ imEnt , imSal ] )
Figure 6.1: Resultado obtenido
12