Smoothing Problem in 2D Images with Tissue-like P Systems and Parallel Implementation

Smoothing Problem in 2D Images with Tissue-like P Systems and Parallel Implementation Francisco Pe˜na-Cantillana1, Daniel D´ıaz-Pernil2, Hepzibah A. Christinal2,3, Miguel A. Guti´errez-Naranjo1 1Research Group on Natural Computing Department of Computer Science and Artificial Intelligence University of Sevilla, Spain [email protected],[email protected] 2Computational Algebraic Topology and Applied Mathematics Research Group Department of Applied Mathematics I University of Sevilla, Spain [email protected] 3Karunya University, Coimbatore, Tamilnadu, India [email protected] Summary. Smoothing is often used in Digital Imagery to reduce noise within an image. In this paper we present a Membrane Computing algorithm for smoothing 2D images in the framework of tissue-like P systems. The algorithm has been implemented by using a novel device architecture called CUDATM , (Compute Unified Device Architecture). We present some examples, compare the obtained time and present some research lines for the future. 1 Introduction The study of digital images [22] has seen a large progress over the last decades. The aim of dealing with an image in its digital form is improving its quality, in some sense, or to simply achieving some artistic effect. The physical properties of camera technology are inherently linked to different sources of noise, so the application of smoothing algorithm are necessary for an appropriate use of the images. Smoothing is often used to reduce such noise within an image. In this paper we use Membrane Computing4techniques for smoothing 2D images in the framework of tissue-like P systems. The algorithm has been implemented by using a novel device architecture called CUDATM , (Compute Unified Device Architecture) [16, 23]. CUDATM is a general purpose parallel computing architecture that allows the parallel NVIDIA Graphics Processors Units (GPUs) 4We refer to [19] for basic information in this area, to [20] for a comprehensive presentation and the web site [24] for the up-to-date information. 318 F. Pe˜na-Cantillana et al. to solve many complex computational problems5in a more efficient way than on a CPU . This parallel architecture has been previously used in Membrane Computing [1, 2, 3] but, to the best of our knowledge, this is the first time that it is used for implementing smoothing algorithms with tissue-like P systems. Dealing with Digital Imagery has several features which make it suitable for techniques inspired by nature. One of them is that it can be parallelized and locally solved. Regardless how large the picture is, the segmentation process can be performed in parallel in different local areas. Another interesting feature is that the basic necessary information can be easily encoded by bio-inspired representations. In the literature, one can find several attempts for bridging problems from Digital Imagery with Natural Computing as the works by K.G. Subramanian et al. [4, 5] or the work by Chao and Nakayama where Natural Computing and Algebraic Topology are linked by using Neural Networks [6] (extended Kohonen mapping). Recently, new approaches have been presented in the framework of Membrane Computing [7, 11]. In [8, 9, 10], Christinal et al. started a new bio-inspired research line where the power and efficiency of tissue-like P systems [12, 13] were applied to topological processes for 2D and 3D digital images. The paper is organised as follows: Firstly, we recall some basics of tissue-like P systems and the foundations of Digital Imagery. Next we present our P systems family and an easy example showing different results by using different thresholds. In Section 3 we present the implementation in CUDATM of the algorithm and show an illustrative example, including a comparative of the time obtained in the different variants. The paper finishes with some final remarks and hints for future work. 2 Preliminaries In this section we provide some basics on the used P system model, tissue-like P systems, and on the foundation of Digital Imagery. Tissue-like P systems [14, 15] have two biological inspirations: intercellular communication and cooperation between neurons. The common mathematical model of these two mechanisms is a network of processors dealing with symbols and communicating these symbols along channels specified in advance. Formally, a tissue-like P system with input of degree q≥1 is a tuple Π= (Γ, Σ, E, w1, . . . , wq,R, iΠ, oΠ), where 1. Γis a finite alphabet, whose symbols will be called objects; 2. Σ(⊂Γ) is the input alphabet; 3. E ⊆ Γis the alphabet of objects in the environment; 4. w1, . . . , wqare strings over Γrepresenting the multisets of objects associated with the cells at the initial configuration; 5For a good overview, the reader can refer to [17, 18]. Smoothing Problem in 2D Images with Tissue-like P Systems 319 5. Ris a finite set of communication rules of the following form: (i, u/v, j) for i, j ∈ {0,1,2, . . . , q}, i 6=j,u, v ∈Γ∗; 6. iΠ∈ {1,2, . . . , q}is the input cell; 7. oΠ∈ {0,1,2, . . . , q}is the output cells A tissue-like P system of degree q≥1 can be seen as a set of qcells (each one consisting of an elementary membrane) labelled by 1,2, . . . , q. We will use 0 to refer to the label of the environment, iΠdenotes the input region and oΠdenotes the output region (which can be the region inside a cell or the environment). The strings w1, . . . , wqdescribe the multisets of objects placed in the qcells of the P system. We interpret that E ⊆ Γis the set of objects placed in the environment, each one of them available in an arbitrary large amount of copies. The communication rule (i, u/v, j) can be applied over two cells labelled by i and jsuch that uis contained in cell iand vis contained in cell j. The application of this rule means that the objects of the multisets represented by uand vare interchanged between the two cells. Note that if either i= 0 or j= 0 then the objects are interchanged between a cell and the environment. Rules are used as usual in the framework of membrane computing, that is, in a maximally parallel way (a universal clock is considered). In one step, each object in a membrane can only be used for one rule (non-deterministically chosen when there are several possibilities), but any object which can participate in a rule of any form must do it, i.e., in each step we apply a maximal set of rules. Aconfiguration is an instantaneous description of the P system Π. Given a configuration, we can perform a computation step and obtain a new configuration by applying the rules in a parallel manner as it is shown above. A sequence of computation steps is called a computation. A configuration is halting when no rules can be applied to it. Then, a computation halts when the P system reaches a halting configuration. Next we recall some basics on Digital Imagery6. Apoint set is simply a topological space consisting of a collection of objects called points and a topology which provides for such notions as nearness of two points, the connectivity of a subset of the point set, the neighborhood of a point, boundary points, and curves and arcs. For a point set Xin Z, a neighborhood function from Xin Z, is a function N:X→2Z. For each point x∈X,N(x)⊆Z. The set N(x) is called a neighborhood for x. There are two neighborhood function on subsets of Z2which are of particular importance in image processing, the von Neumann neighborhood and the Moore neighborhood. The first one N:X→2Z2is defined by N(x) = {y:y= (x1±j, x2)or y = (x1, x2±k), j, k ∈ {0,1}}, where x= (x1, x2)∈X⊂Z2. While the Moore neighborhood M:X→2Z2is defined by M(x) = {y:y= (x1±j, x2±k), j, k ∈ {0,1}}, where x= (x1, x2)∈X⊂Z2. The von Neumann 6We refer the interested reader to [21] for a detailed introduction. 320 F. Pe˜na-Cantillana et al. and Moore neighborhood are also called the four neighborhood (4-adjacency) and eight neighborhood (8-adjacency), respectively. An Z-valued image on Xis any element of ZX. Given an Z-valued image I∈ZX, i.e. I:X→Z, then Zis called the set of possible range values of I and Xthe spatial domain of I. The graph of an image is also referred to as the data structure representation of the image. Given the data structure representation I={(x, I(x)) : x∈X}, then an element (x, I(x)) is called a picture element or pixel. The first coordinate xof a pixel is called the pixel location or image point, and the second coordinate I(x) is called the pixel value of Iat location x. For example, Xcould be a subset of Z2where x= (i, j) denotes spatial location, and Zcould be a subset of Nor N3, etc. So, given an image I∈ZZ2, a pixel of Iis the form ((i, j), I(x)), which be denoted by I(x)ij. We call the set of colors or alphabet of colors of I,CI⊆Z, to the image set of the function Iwith domain Xand the image point of each pixel is called associated color. We can consider an order in this set. Usually, we consider in digital image a predefined alphabet of colors C ⊆ Z. We define h=|C| as the size (number of colors) of C. In this paper, we work with images in grey scale, then C={0, . . . , 255}, where 0 codify the black color and 255 the white color. Aregion could be defined by a subset of the domain of Iwhose points are all mapped to the same (or similar) pixel value by I. So, we can consider the region Rias the set {x∈X:I(x) = i}but we prefer to consider a region ras a maximal connected subset of a set like Ri. We say two regions r1, r2are adjacent when at less a pair of pixel x1∈r1and x2∈r2are adjacent. We say x1and x2are border pixels. If I(x1)< I(x2) we say x1is an edge pixel. The set of connected edge pixels with the same pixel value is called a boundary between two regions. The purpose of image enhancement is to improve the visual appearance of an image, or to transform an image into a form that is better suited for human interpretation or machine analysis. There exists a multitude of image enhancement techniques, as are averaging of multiple images, local averaging, Gaussian smoothing, max-min sharpening transform, etc. One of the form to enhancement an image could be eliminate non important regions of an image, i.e., remove regions which do not provide relevant information. This technique is known as smoothing. 2.1 A Family of Tissue-like P Systems Given a digital image with n2pixels and n∈Nwe define a tissue-like P system whose input is given by the pixels of the image encoded by the objects aij , where 1≤i, j ≤nand a∈ C. Next, we shall give some outlines how to prove that our smoothing problem can be solved in a logarithmic number of steps using a family of tissue-like P systems Π. We define a family of tissue-like P systems to do an smoothing of a 2D image. For each image of size n2with n∈N, we consider the tissue-like P system with input of degree 1: Smoothing Problem in 2D Images with Tissue-like P Systems 321 Π(r, n) = (Γ, Σ, E, w1,R, iΠ, oΠ), where •Γ=Σ∪ E, •Σ={aij :a∈ C,1≤i, j ≤n}, • E ={aij : 1 ≤i, j ≤n, a ∈ C}, •w1=∅, •Ris the following set of communication rules: – (1, aijbkl/aijakl,0), ·for 1 ≤i, j ≤n, ·a, b ∈ C,a < b and d(a, b)≤r, These rules are used to simplify the image. If we have two colors whose distance in the alphabet of colors of the image is small, then we change the high color by the small one. Of this manner, we change the regions structure. •iΠ=oΠ= 1. Each P system works as follows: We take pairs of adjacent pixels and change the color of the pixel with lower color. We do it in a parallel way with all the possible pairs of pixels. In the next step, we will repeat the previous process, but the colors of the pixels have could to be changed. So, we need, in the worst case, a linear number of steps to do all the possible changes to obtain an smoothing of our image. Fig. 1. An example 2.2 An Easy Example In this section, we show the results obtained by the application of our method. Our input image (of size 30×30) can be seen in Figure 1. In this case, 0 represents 322 F. Pe˜na-Cantillana et al. to white color and 30 represents to black color, i.e., the inverse order for Cis considered. Fig. 2. Theoretical Smoothness of the Figure 1. Working with different thresholds provides different results as we can observe in Figure 2. If we take r= 5, then we get the first image, and when the threshold is r= 10 the second image is obtained. By using this method, it is clear that the structure of regions is (more or less) conserved with a threshold of r= 5 and we need to a high number to obtain a more simplified image. Moreover, we can observe in both cases the color of regions is similar to the regions of input image in this method. Therefore, we can use this technique for smoothing images, and clarify the structure of a region without eliminate important information. Bearing in mind the size of the input data is O(n2), |C| =hand ris the threshold used with both membrane solutions. The amount of necessary resources for defining the systems of our families and the complexity of our problem is determined in the following table: Smoothing Problem Complexity Dynamical Number of steps of a computation O(n) Necessary Resources Size of the alphabet n2·h Initial number of cells 1 Initial number of objects 0 Number of rules O(n2·h) Upper bound for the length of the rules 4 3 Parallel Implementation GPUs constitute nowadays a solid alternative for high performance computing, and the advent of CUDATM allow programmers a friendly model to accelerate a broad range of applications. The way GPUs exploit parallelism differ from multi-core Smoothing Problem in 2D Images with Tissue-like P Systems 323 CPUs, which raises new challenges to take advantage of its tremendous computing power. GPU is especially well-suited to address problems that can be expressed as data-parallel computations. GPUs can support several thousand of concurrent threads providing a massively parallel environment. This parallel computation model leads us to look for a highly parallel computational technology where a parallel simulator can run efficiently. In this paper, we present a parallel software tool based in our membrane solution for smoothing images. It has been developed by using Microsoft Visual Studio 2008 Professional Edition (C++) with the plugging Parallel Nsight (CUDATM ) under Microsoft Windows 7 Professional with 32 bits. To implement the P systems, CUDATM C, an extension of C for implementations of executable kernels in parallel with graphical cards NVIDIA has been used. It has been necessary the nvcc compiler of CUDATM Toolkit. Moreover, we use libraries from openCV to the treatment of input and output images. Microsoft Visual Studio 2008 is responsible for calling to the compilers to build the objects, and to link them with the final program. This allows us to deal with images stored in .BMP, .DIB, .JPEG, .JPG, .JPE, .PNG, .PBM, .PGM, .PPM, .SR, .RAS, .TIFF and .TIF formats The experiments have been performed on a computer with a CPU Intel Pentium 4 650, with support for HT technology which allows to work like two CPUs of 32 bits to 3412 MHz. Computer has 2 MB of L2 cache memory and 1 GB DDR SDRAM of main memory with 64 bits bus wide to 200 MHz. Moreover, it has a hard disc of 160 GB SATA2 with a transfer rate of 300 Mbps in a 8 MB buffer. The graphical card (GPU) is an NVIDIA Geforce 8600 GT composed by 4 Stream Processors with a total of 32 cores to 1300 MHz and executes 512 threads per block as maximum. It has a 512 MB DDR2 main memory, but 499 MB could be used by processing in a 128 bits bus to 700 MHz. So, the transfer rate obtained is by 22.4 Gbps. For constant memory used 64 KB and for shared memory 16 KB (It is not a good data for a good CUDATM graphical card). Its Compute Capability is 1.1 (from 1.0 to 2.1), then we can obtain a lot of improvements in the efficiency of the algorithms. If we compare CPU and GPU, we observe that the former has two cores to 3412 MHz and the latter has 32 to 1300Mhz and has larger memory. We have developed two applications of our P systems. In this case, we consider the natural order in the set of colors C={0, . . . , 255}. In the first one, we have considered a deterministic implementation, where we work with the Moore neighborhood. So, the system checks if the rules can be applied for eight adjacent pixels. In the second one, we have considered an random selection of an adjacent pixel to work. Then, the system checks only one possibility chosen in an random way. Of this form, we simulate the characteristic non determinism of P systems with random. Moreover, we have decided to stop the system before the halting configuration, because more than an appropriate number of parallel steps of processing could be non-operative. In fact, in the deterministic version, the process could fin- 324 F. Pe˜na-Cantillana et al. ish before the pre-fixed number of steps. So, the system needs some time to check this possibility. In the second one, it is not necessary to look at this question. Fig. 3. Original image Fig. 4. Deterministic version, Threshold 50: 0, 5, 10, 20, 32 and 44 steps, respectively. Fig. 5. Deterministic version. First: Threshold 75, step 193. Second: Threshold 125, step 193. Smoothing Problem in 2D Images with Tissue-like P Systems 325 We consider the image of size 640 ×400 in Figure 3. When we take the deterministic application of our software, we can check that if we use an threshold r= 50, our software smooths the original image (see Figure 4) using 44 parallel steps. Nonetheless, when we work with a higher threshold, new important regions are changed, and the output image is different. (See the images of Figure 5). When we consider the random version of our software, we can check that if we use an threshold r= 50 we need 300 steps, but the differences with the resulting image with 150 or 200 steps are minimum, as we can see in Figure 6. When we take higher thresholds, as in the Figure 7, we can check that new regions change of colors. Fig. 6. Random version, Threshold 50: 0, 50, 100, 150, 200 and 300 steps, respectively. Fig. 7. Random version.First: Threshold 75, step 1000. Second: Threshold 125, step 800. We present a study the running time of our software for both cases with different thresholds, showed in the above examples, with the next table. We can observe that deterministic version of our software needs less time with respect to the random version. In the first one, it applies eight rules for each pixels while, in the second one, it applies only one rule for each pixel. Moreover, we need a running time to implement the random in each step for each pixel.