Improvement of Modal Matching Image Objects in Dynamic Pedobarography using Optimization Techniques

Electronic Letters on Computer Vision and Image Analysis 5(3):1-20, 2005 Improvement of Modal Matching Image Objects in Dynamic Pedobarography using Optimization Techniques João Manuel R. S. Tavares* and Luísa Ferreira Bastos+ * Laboratório de Óptica e Mecânica Experimental, Instituto de Engenharia Mecânica e Gestão Industrial / Departamento de Engenharia Mecânica e Gestão Industrial, Faculdade de Engenharia da Universidade do Porto, Rua Dr. Roberto Frias, s/n, 4200-465, Porto, PORTUGAL + Laboratório de Óptica e Mecânica Experimental, Instituto de Engenharia Mecânica e Gestão Industrial, Rua Dr. Roberto Frias, s/n, 4200-465, Porto, PORTUGAL Received 20 December 2004; accepted 9 March 2005 Abstract This paper presents an improved approach for matching objects represented in dynamic pedobarography image sequences, based on finite element modeling, modal analysis and optimization techniques. In this work, the determination of correspondences between objects data points is improved by using optimization techniques and, because the number of data points of each object is not necessary the same, a new algorithm to match the excess points is also proposed. This new matching algorithm uses a neighbourhood criterion and can overcome some disadvantages of the usual “one to one” matching. The considered approach allows the determination of correspondences between 2D or 3D objects data points, and is here apply in dynamic pedobarography images. Key Words: Computacional Vision, Deformable Objects, Dynamic Pedobarography, FEM, Matching, Medical Imaging, Modal Analysis, Optimization Techniques. 1 Introduction In several areas of Computational Vision, one of the main problems consists in the determination of correspondences between objects represented in different images, and on the computation of robust canonical descriptors that can be used for their recognition. In this paper, is presented a methodology to address the above problem, based on an approach initially proposed by Sclaroff [1], [2], and in this work improved by using optimization algorithms in the matching phase. A new algorithm to determine the correspondences between excess models nodes is also proposed, and can be used when the objects to be matched are represented by different number of data points. With this new algorithm, we can successfully overcome some disadvantages of the usual “one to one” matching as, for example, loss of information along image sequences. The application of the proposed matching methodology on deformable objects represented in dynamic pedobarography image sequences leads to very promising results, and will be discussed in this paper. Correspondence to: <[email protected]> Recommended for acceptance by Perales F., Draper B. ELCVIA ISSN: 1577-5097 Published by Computer Vision Center / Universitat Autonoma de Barcelona, Barcelona, Spain 2 Tavares et al. / Electronic Letters on Computer Vision and Image Analysis 5(3):1-20, 2005 The following sections present a brief introduction to the background problem, the dynamic pedobarography principle, the used object models, the proposed matching methodology, experimental results obtained on deformable objects represented in dynamic pedobarography images and some conclusions. 1.1 Background There is an eigen methods class [2] that derives its parameterization directly from objects data shape. Some of these techniques also try to determine, explicitly and automatically, the correspondences between characteristic points’ sets, while others try to match images using more global approaches instead of local ones. Each eigen method decomposes the object deformation in an orthogonal and ordered base. Usually, solution methods for the matching problem include several restrictions that prevent inadequate matches according to some criteria, as for example: order [3], [4]; rigidity [3], [4]; unicity [5]; visibility [6]; and proximity [2]. Some of these methods are image correlation (it is presumed that the images are similar) [7], point proximity [8], and smoothness of disparity fields [3]. The matching problem can also be interpreted as an optimization problem, in which the objective function can, for example, depend on any criteria mentioned in the previous paragraph, and the restrictions considered must form a non-empty space of possible solutions. To solve this optimization problem, it can be used dynamic programming [3], graphs [4] and convex minimization [7]. Non-optimal approaches include, for example, greedy algorithms [9], simulated annealing [10], relaxation [5], etc. To determine correspondences between two objects, Belongie [11] considered shapes context and a similar optimization technique to the one used in this work. Although shape description algorithms have usually a higher computational efficiency, the modeling methodology considered in this work as the major advantage of attributing a physical behaviour to each object to be matched, through the consideration of a virtual elastic material. 2 Dynamic Pedobarography Pedobarography refers to measuring and visualizing the distribution of pressure under the foot sole. The recording of pedobarographic data along the duration of a step, in normal walking conditions, permits the dynamic analysis of the foot behavior. This introduction of the time dimension augments the potential of this type of clinical examination as an auxiliary tool for diagnostics and therapy planning [12]. The basic pedobarography system consists of a transparent plate trans-illuminated through its borders in such a way that the light is internally reflected. The plate is covered on its top by a single or dual thin layer of soft porous plastic material where the pressure is applied (see Figure 1). reflected light glass or acrylic plate foot pressure opaque layer lamp foot pressure Figure 1: Basic (pedo)barography principle. When observed from below, in the absence of applied pressure the plate is dark. However, when pressure is applied on top of the plastic layer, the plate displays bright areas that correspond to the light crossing the plate after reflection on the plastic layer. This reflection occurs due to the alteration of the local relation of light refraction indices resulting from the depletion of the air interface between the glass plate and the plastic layer. A good choice of materials and an adequate calibration of the image acquisition system, allow a nearly proportional relation between the local pressure and the observed brightness. Using a practical setup as the one represented in Figure 2, a time sequence of pressure images can be captured. Figure 3 shows thirteen images of a captured sample sequence; as can be verified, the image data is Tavares et al. / Electronic Letters on Computer Vision and Image Analysis 5(3):1-20, 2005 3 very dense, as opposed to other measuring methods, and very rich in terms of the information it conveys on the interaction between the foot sole and the flat plate. Video 1 shows the image sequence considered. pedobarography table camera mirror glass + contact layer Figure 2: Basic setup of a pedobarography system. 0 1 2 3 4 5 6 7 8 9 10 11 12 Figure 3: Example of a captured sequence composed by thirteen images. 4 Tavares et al. / Electronic Letters on Computer Vision and Image Analysis 5(3):1-20, 2005 Video 1: Animation of the thirteen images of the sample sequence (click to see the video). 3 Object Models In the initial stages of the work, the object contours in each image were extracted and the matching process was oriented to the contours’ pixels [2], [13]. A practical difficulty arising from this approach is the possible existence of more than one contour for the object represented in each image (i. e. see Figure 7). To find the correspondence between each contours pair along the image sequence, two possible solutions were considered: i) use of a Kalman filtering (see [13], for example) approach to estimate and track the location of the contours’ centroids along the image sequence; ii) use of a measure of the deformation energy necessary to align each contour pairs, selecting the lower energy pairs. However, an additional problem is still present: the possibility that along the image sequence various contours will merge or split. In order to accommodate this possibility, a new model has been developed, similar to the one used in various applications working with controlled environment, such as in face analysis and recognition [14], [15]: The brightness level of each pixel is considered as the third coordinate of a 3D surface point. The resulting single surface model solves the two aforementioned problems. The use of the surface model, also simplifies the consideration of isobaric contours, which are important in pedobarographic analysis, either for matching contours of equal pressure along the time sequence or for matching contours of different pressure in a single image [12], [13]. The following sections describe the object models used in this work and briefly describe their construction. Each model has its own advantages and shortcomings; for every particular problem, the best choice must be made [12], [13]. 3.1 Contour Model To determine the correspondence among two contours were used two modeling approaches: • For each contour is used a single 2D Sclaroff’s isoparametric finite element. In building this type of element, no previous ordering of the data points is required and Gaussian functions are used as interpolations functions. The method to determine the mass and stiffness matrices for this 2D element is described in, for example, [1]. • For each contour is built a finite elements model using linear axial finite elements (Figure 4). For this type of discretisation a previous ordering of the contour data points is required. The matrix formulation for these finite elements can be found in [16], for example. Tavares et al. / Electronic Letters on Computer Vision and Image Analysis 5(3):1-20, 2005 5 3 e1e2 e3 e4 e5 e6 e7 e81 2 4 5 6 7 8 Figure 4: Modeling a contour by using a set of axial finite elements. i e To determine in each image the contour pixels, are used standard image processing and analysis techniques; namely, thresholding, edge enhancement, hysteresis line detection and tracking [2]. For example, Figures 6 and 7 show the intermediate result and the final contours obtained from the image of Figure 5. The sampled contours obtained from the thirteen images shown in Figure 3 are presented in Figure 8. Figure 5: Image (negated) where contours must be found. Figure 6: Result image after edge enhancement. Figure 7: Contours obtained by a line detection and tracking algorithm [2]. 3.2 Surface Model For the surface model, two approaches were also considered: • A single 3D Sclaroff’s isoparametric finite element is used for each surface to be matched. Again, it must be noticed that there is no requirement for previous ordering of the surface nodes. The matrix building for these finite elements can be found in [1]. • Each surface model is built by using linear axial finite elements (Figure 9). The previous ordering of the surface data points is required. The matrix formulation for these finite elements can be found in [16], for example. The used methodology to determine the data points of each surface can be summarized as follows: 1. Noise pixels (that is, pixels with brightness lower than a calibration threshold level) are removed and a Gaussian-shaped smoothing filter is applied to the image (Figure 10); 2. The circumscribing rectangle of the object to be modeled is determined and the image is sampled within that area (Figure 11); 3. A 2D Delaunay triangulation (see, for example, [25, 26]) is performed on the sampled points, using the point brightness as the third coordinate; 4. In order to reduce the number of nodes used, and thus the computational cost, the triangular mesh is simplified using a decimation algorithm (see, for example, [25, 26]); 5. To reduce the high frequency noise associated to the mesh, is used a Laplacian smoothing algorithm (see, for example, [25, 26]); 6. Finally, in order to have similar ranges of values in all coordinates, a scale change is performed on the third coordinate (derived from brightness) (Figure 12). 6 Tavares et al. / Electronic Letters on Computer Vision and Image Analysis 5(3):1-20, 2005 0 1 2 3 4 5 6 7 8 9 10 11 12 Figure 8: Sampled contours obtained from the original image sequence of Figure 3. Figure 9: Modeling a surface by using a set of axial finite elements. (Each node is connected to its neighbors through axial elements.) Tavares et al. / Electronic Letters on Computer Vision and Image Analysis 5(3):1-20, 2005 7 Figure 10: Image (negated) after noise removal and Gaussian filtering. Figure 11: Object sampling. Figure 12: Resulting surface. The surfaces obtained from the original images presented in Figure 3 are visible in Figure 13. The original images with identification (ID) 0 (zero) and 1 (one) were not considered. 2 3 4 5 6 7 8 9 10 11 12 Figure 13: Surfaces obtained from the last eleven images of the example sequence. 3.3 Isobaric Contour Model As in the two previous models, the approaches used to match isobaric contours are: 8 Tavares et al. / Electronic Letters on Computer Vision and Image Analysis 5(3):1-20, 2005 • A single Sclaroff’s isoparametric finite element, either 2D or 3D, is used to model each contour. • A set of linear axial finite elements are used to build each contour model. The isobaric contours are extracted from the correspondent surface, beforehand obtained using the procedures described in the previous section (Figure 14). Figure 14: Ten isobaric contours extracted from the surface of figure 12. 4 Matching Methodology Figure 15 displays a diagram of the adopted physical matching methodology. The locations of the objects data points in each image, [ ] 1 =…m XXX , are used as the nodes of a finite elements model made of an elastic material. Next, the eigenmodes { } i φ of the model are computed, providing an orthogonal description of the object and its natural deformations, ordered by frequency. Using a matrix based notation, the eigenvectors matrix [ ] Φ and the eigenvalues diagonal matrix [ ] Ω can be written as in equation (1) for 2D objects and as in equation (2) for 3D objects. The eigenvectors, also called shape vectors [1], [2], [16], [17], describe how each vibration mode deforms the object by changing the original data point locations: { } =+ deformed i XXa φ . [] {} {} {} {} {} {} [] 1 2 1 12 2 12 0 and 0 T T m T m m T m u u v v ω φφ ω ⎡⎤ ⎢⎥ ⎢⎥ ⎡ ⎤ ⎢⎥ ⎢ ⎥ ⎢⎥ ⎡⎤ Φ= = Ω= ⎢ ⎥ ⎣⎦ ⎢⎥ ⎢ ⎥ ⎢⎥ ⎣ ⎦ ⎢⎥ ⎢⎥ ⎢⎥ ⎣⎦     , (1) [] {} {} {} {} {} {} {} {} [] 1 2 1 1 13 2 3 1 0 and 0 T T m T m T m m T T m u u v v w w ω φφ ω ⎡⎤ ⎢⎥ ⎢⎥ ⎢⎥ ⎢⎥ ⎢⎥ ⎡ ⎤ ⎢⎥ ⎢ ⎥ ⎢⎥ ⎡⎤ Φ= = Ω= ⎢ ⎥ ⎣⎦ ⎢⎥ ⎢ ⎥ ⎢⎥ ⎣ ⎦ ⎢⎥ ⎢⎥ ⎢⎥ ⎢⎥ ⎢⎥ ⎢⎥ ⎣⎦      . (2) The first three (in 2D) or six (in 3D) vibration modes are the rigid body modes of translation and rotation; the remaining modes are non-rigid [1], [2], [16], [17]. In general, lower frequency modes describe global deformations, while higher frequency modes essentially describe local deformations. This type of ordering, from global to local behaviour, is quite useful for object matching and recognition. Tavares et al. / Electronic Letters on Computer Vision and Image Analysis 5(3):1-20, 2005 9 Figure 15: Diagram of the adopted matching methodology. The eigenmodes also form an orthogonal, object-centred coordinate system for the location of the data points, i.e., the location of each point is uniquely described in terms of each eigenmode displacement. The transformation between the Cartesian image coordinates and the modal coordinates system is achieved through the eigenvectors matrix of the physical model. Two sets of data points, for example, corresponding to the objects represented in two different images of a sequence, are to be compared in the modal eigenspace. The main idea is that the low order modes of two similar objects will be very close, even in the presence of an affine transformation, a non-rigid deformation, a local shape variation, or noise. Using the above concept, data correspondence is obtained by modal matching. In this work, two matching search procedures are considered: (i) a local search or (ii) a global search. Both search strategies consider an affinity matrix [ ] Z , constructed from the Euclidian distance between the characteristic vectors of each physical model, whose elements are, for 2D and for 3D, respectively: {} {} {} {} 2 ,1, ,1,+ =− +− ij ti tj ti tj Zu u v v2 + , (3) {} {} {} {} {} {} 22 ,1, ,1, , 1,++ =− +− +− ij ti t j ti t j ti t j Zu u v v w w2 + . (4) The local search strategy was proposed in the original modal matching methodology [1], [2], [13], and basically it consists in seeking each row of the affinity matrix for its lowest value, considering the associated correspondence if, and only if, that value is also the lowest of the related column. With this search philosophy, the models nodes are considered as independent entities and so the original object structure is ignore. 16 Tavares et al. / Electronic Letters on Computer Vision and Image Analysis 5(3):1-20, 2005 Local matching: 52 Optimization matching: 76 Figure 23: Matches obtained between isocontours with ID 1 (76 nodes) and 2 (76 nodes). Local matching: 40 Optimization matching: 70 Excess nodes matching: 76 Figure 24: Matches obtained between isocontours with ID 4 (76 nodes) and 5 (70 nodes). 6 Conclusions The several experimental tests carried through this work, some reported in this paper, allow the presentation of some observations and conclusions. The physical methodology proposed, for the determination of correspondences between two objects, using optimization techniques on the matching phase, when compared with the one previously developed that considers local search, obtained always an equal or higher number of satisfactory matches. It was also verified that the number of matches found is independent from the optimization algorithm considered. In some experimental cases, in order to obtain a higher number of satisfactory matches using the local search strategy, the parameters of the physical methodology had to be carefully chosen. In those same cases, the application of the optimization strategy in the matching phase, beyond the good quality of the matches found, revealed less sensible to the methodology parameters. This suggests that using the proposed global Tavares et al. / Electronic Letters on Computer Vision and Image Analysis 5(3):1-20, 2005 17 search strategy in the matching phase, the adopted physical methodology become easier to handle and also more adjustable to different kinds of applications. Local matching: 24 Optimization matching: 46 Excess nodes matching: 54 Figure 25: Matches obtained between isocontours with ID 6 (54 nodes) and 7 (46 nodes). Local matching: Optimization matching: Excess nodes matching: Figure 26: Matches obtained between the eleven isocontours extracted from the image with ID 6 of the example sequence. 18 Tavares et al. / Electronic Letters on Computer Vision and Image Analysis 5(3):1-20, 2005 Local Optimization Excess Nodes Isocontours ID Nº Nodes Nº M atch % Match Nº Match atch % Match Nº M 1, 2 76/76 52 76 68,4% 100% 2, 3 76 4 58 ,4% /7 78 74 100% 76 3, 4 74/76 50 67,6% 74 100% 76 4, 5 76/70 40 57,1% 70 100% 76 5, 6 70/54 32 59,3% 54 100% 70 6, 7 54/46 24 52,2% 46 100% 54 7, 8 46/38 22 57,9% 38 100% 46 8, 9 38/34 17 50,0% 34 100% 38 9, 10 34/34 32 94,1% 34 100% 1 0, 11 34/34 17 50,0% 34 100% Table mber (N h) and centag ches (% Match) contou ted fr the im ID 6 e exa quence. 3: Nu º Matc per e of mat obtained with the eleven iso rs extrac om age with of th mple se Local matching Optimization matching Videos 5 and 6: Matches found between the isocontours extracted from the image with ID 6 of the example sequence (click to see the video). To ocal search strategy, the global search r of eigenvectors in the affinity matrix construction. This suggests that the total computational effort of the global matching methodology can rithms considered. The execution time of LAPm alg nclude that the referred algorithm can become an interesting base for the development of new sol have satisfactory matching results in some of the examples considered, when compared with the l strategy always required an inferior numbe be reduced if optimization techniques are considered. In terms of execution time, the optimization algorithm that uses the Hungarian method showed to be of low efficiency. In the several experimental examples performed, the Simplex algorithm for flow problems revealed the most efficient among the optimization algo orithm was higher than the Simplex algorithm for flow problems, even being a more specific algorithm for this type of problems. This can be dew to the interval in which lay the elements of the affinity matrix, [0; 1], since when this algorithm was tested in [27] it revealed the most efficient when the considered interval was [1; 100]. In the several experimental tests performed considering contour objects, the implemented algorithm for the determination of correspondences of the excess data points always finds satisfactory matches. That allows us to co utions for the determination of matches of type “one to many” and vice versa, and that should be ported to more complex objects (i. e., surfaces). The experimental results shown in this paper, confirm that the proposed physical methodology can satisfactory match objects represented in dynamic pedobarography images, and that the use of the pixel Tavares et al. / Electronic Letters on Computer Vision and Image Analysis 5(3):1-20, 2005 19 bri [1] ff and A. Pentland, “Modal Matching for Correspondence and Recognition”, IEEE on Pattern Analysis and Machine Intelligence, vol. 17, pp. 545-561, 1995 dia, 1998 ca [7] attern Analysis and Machine Intelligence, vol. 25, pp. 187-199, 2003 nition Letters, vol. 16, pp. 23-31, 1995 nd Machine Intelligence, vol. 24, pp. 509-522, 2002 (RecPad'00), [14] chine Vision Conference (BMVC’94), 1994 1, 1996 [17] 93 espondence of Objects ering (BioEng'2003), [19] Espaço Modal”, in Faculdades de Engenharia e Ciências: Universidade do Porto, [21] hnik Berlin, Division Scientific Computing, Department Optimization, 2000 ghtness values as a third Cartesian coordinate is very satisfactory, both in terms of its interpretation as pressure, and in solving the problems associated to merging or the splitting of objects. References S. E. Sclaro Transactions [2] J. M. R. S. Tavares, PhD Thesis, “Análise de Movimento de Corpos Deformáveis usando Visão Computacional”, in Faculdade de Engenharia da Universidade do Porto, Portugal, 2000 [3] Y. Ohta and T. Kanade, “Stereo by Intraand Inter-Scanline Search using Dynamic Programming”, IEEE Transactions on Pattern Analysis and Machine Intelligence, vol. 7, pp. 139-154, 1985 [4] S. Roy and I. J. Cox, “A Maximum-Flow Formulation of the N-camera Stereo Correspondence Problem”, presented at International Conference on Computer Vision (ICCV'98), Bombay, In [5] S. Gold, A. Rangarajan, C. P. La, S. Pappu and E. Mjolsness, “New algorithms for 2D and 3D point matching: pose estimation and correspondence”, Pattern Recognition, vol. 31, pp. 1019-1031, 1998 [6] C. Silva, PhD Thesis, “3D Motion and Dense Structure Estimation: Representation for Visual Perception and the Interpretation of Occlusions”, in Instituto Superior Técnico: Universidade Técni de Lisboa, 2001 J. L. Maciel and J. P. Costeira, “A Global Solution to Sparse Correspondence Problems”, IEEE Transactions on P [8] Z. Zhang, “Iterative Point Matching for Registration of Free-Form Curves”, INRIA, Technical Report RR-1658, April 1992 [9] M. S. Wu and J. J. Leou, “A Bipartite Matching Approach to Feature Correspondence in Stereo Vision”, Pattern Recog [10] J. P. P. Starink and E. Backer, “Finding Point Correspondences using Simulated Annealing”, Pattern Recognition, vol. 28, pp. 231-240, 1995 [11] S. Belongie, J. Malik and J. Puzicha, “Shape Matching and Object Recognition using Shape Context”, IEEE Transactions on Pattern Analysis a [12] A. J. Padilha, L. A. Serra, S. A. Pires and A. F. N. Silva, “Caracterização Espacio-Temporal de Pressões Plantares em Pedobarografia Dinâmica”, FEUP/INEB, Relatório Interno, 1995 [13] J. M. R. S. Tavares, J. Barbosa and A. Padilha, “Matching Image Objects in Dynamic Pedobarography”, presented at 11th Portuguese Conference on Pattern Recognition Porto, Portugal, 2000 T. F. Cootes and C. J. Taylor, “Modelling Object Appearance Using The Grey-Level Surface”, presented at British Ma [15] B. Moghaddam, C. Nastar and A. P. Pentland, “Bayesian Face Recognition using Deformable Intensity Surfaces”, MIT Media Laboratory, Technical Report Nº 37 [16] K.-J. Bathe, Finite Element Procedures, Prentice Hall, 1996 S. Graham Kelly, Fundamentals of Mechanical Vibrations, McGraw-Hill, 19 [18] L. F. Bastos and J. M. R. S. Tavares, “Optimization in Modal Matching for Corr Nodal Points”, presented at 7th Portuguese Conference on Biomedical Engine Fundação Calouste Gulbenkian, Lisboa, Portugal, 2003 F. S. Hillier and G. J. Lieberman, Introduction to Operations Research. Mcgraw-Hill International Editions, 1995 [20] L. F. Bastos, MSc Thesis, “Optimização da Determinação das Correspondências entre Objectos Deformáveis no Portugal, 2003 A. Löbel, “MFC - A Network Simplex Implementation”, Konrad-Zuse-Zentrum für Informationstec 20 Tavares et al. / Electronic Letters on Computer Vision and Image Analysis 5(3):1-20, 2005 [22] A. Volgenant, “Linear and Semi-Assignment Problems: A Core Oriented Approach”, Computers and Operations Research, vol. 23, pp. 917-932, 1996 [23] J. M. R. S. Tavares, J. Barbosa and A. Padilha, “Apresentação de um Banco de Desenvolvimento e Ensaio para Objectos Deformáveis”, Revista Electrónica de Sistemas de Informação, vol. 1, 2002 [25] alization Toolkit, 3rd Edition, Kitware Inc, 2002 and Codes for Dense Assignment Problems: The State of [24] R. Davies, Newmat, A matrix library in C++, 2005 The VTK User’s Guide, Kitware Inc., 2003 [26] W. Schroeder, K. Martin and B. Lorensen, The Visu [27] M. Dell’Amico and P. Tooth, “Algorithms The Art”, Discrete Applied Mathematics, 100, pp. 17-48, 2000