High Performance Computing in Remote Sensing - Chapter 5

Chapter 5 Parallel Implementation of the Recursive Approximation of an Unsupervised Hierarchical Segmentation Algorithm James C. Tilton, NASA Goddard Space Flight Center Contents 5.1 Introduction ............................................................ 97 5.2 Description of the Hierarchical Segmentation (HSEG) Algorithm ......... 99 5.3 The Recursive Formulation of HSEG ................................... 100 5.4 The Parallel Implementation of RHSEG ................................ 102 5.5 Processing Time Performance .......................................... 104 5.6 Concluding Remarks ................................................... 106 References .................................................................. 106 The hierarchical image segmentation algorithm (referred to as HSEG) is a hybrid of hierarchical step-wise optimization (HSWO) and constrained spectral clustering that produces a hierarchical set of image segmentations. HSWO is an iterative approach to region growing segmentation in which the optimal image segmentation is found at NR regions, given a segmentation at NR + 1 regions. HSEG’s addition of con-strained spectral clustering makes it a computationally intensive algorithm, for all but the smallest of images. To counteract this, a computationally efficient recursive approximation of HSEG (called RHSEG) has been devised. Further improvements in processing speed are obtained through a parallel implementation of RHSEG. This chapter describes this parallel implementation and demonstrates its computational efficiency on a Landsat Thematic Mapper test scene. 5.1 Introduction Image segmentation is the partitioning of an image into related sections or regions. Forremotelysensedimagesoftheearth,anexampleofanimagesegmentationwould be a labeled map that divides the image into areas covered by distinct earth surface 97 © 2008 by Taylor & Francis Group, LLC 98 High-Performance Computing in Remote Sensing coverssuchaswater,snow,typesofnaturalvegetation,typesofrockformations,types of agricultural crops, and types of other man created development. In unsupervised image segmentation, the labeled map may consist of generic labels such as region 1, region 2, etc., which may be converted to meaningful labels by a post-segmentation analysis. Segmentation is a key first step for a number of approaches to image analysis and compression. In image analysis, the group of image points contained in each region provides a good statistical sampling of image values for more reliable labeling based on region mean feature values. In addition, the region shape can be analyzed for additional clues to the appropriate labeling of the region. In image compression, the regions form a basis for compact representation of the image. The quality of the pre-requisite image segmentation is a key factor in determining the level of performance for these image analysis and compression approaches. Asegmentationhierarchyisasetofseveralimagesegmentationsofthesameimage at different levels of detail in which the segmentations at coarser levels of detail can be produced from simple merges of regions at finer levels of detail. This is useful for applications that require different levels of image segmentation detail depending on theparticularimageobjectssegmented.Auniquefeatureofasegmentationhierarchy that distinguishes it from most other multilevel representations is that the segment or region boundaries are maintained at the full image spatial resolution for all levels of the segmentation hierarchy. In a segmentation hierarchy, an object of interest may be represented by multi-ple image segments in finer levels of detail in the segmentation hierarchy, and may be merged into a surrounding region at coarser levels of detail in the segmentation hierarchy.Ifthesegmentationhierarchyhassufficientresolution,theobjectofinterest will be represented as a single region segment at some intermediate level of segmen-tationdetail.Thesegmentationhierarchymaybeanalyzedtoidentifythehierarchical level at which the object of interest is represented by a single region segment. The object may then be identified through its spectral and spatial characteristics. Addi-tional clues for object identification may be obtained from the behavior of the image segmentations at the hierarchical segmentation levels above and below the level at which the object of interest is represented by a single region. Segmentation hierarchies may be formed through a region growing approach to image segmentation. In region growing, spatially adjacent regions iteratively merge through a specified merge selection process. Hierarchical step-wise optimization (HSWO) is a form of region growing segmentation in which the iterations consist of finding the best segmentation with one region less than the current segmentation [1, 2, 3]. The best segmentation is defined through a mathematical criterion such as a minimum vector norm or minimum mean squared error. An augmentation of HSWO, called HSEG (for hierarchical segmentation), was introduced by Tilton [4] in which an option is provided for merging spatially non-adjacent regions as controlled by a threshold based on previous merges of spatially adjacent regions. This can be thought of as a form of constrained spectral clustering. The introduction of constrained spectral clustering in HSEG makes it a computa-tionally intensive algorithm, for all but the smallest of images. This is because of a © 2008 by Taylor & Francis Group, LLC Parallel Implementation of the Recursive Approximation 99 requirement to evaluate the dissimilarity between all pairs of regions, rather than just spatially adjacent regions. For a 1024 × 1024 pixel image, this leads to the order of 106 dissimilarity evaluations per iteration in the initial processing stages. This computational difficulty is overcome by a recursive approximation of HSEG, called RHSEG. An early version of RHSEG was discussed in Tilton [5]. This re-cursive formulation not only limits the number of comparisons between spatially non-adjacent regions to a more reasonable number, but also lends itself to a straight-forward and efficient implementation on parallel computing platforms. The current parallel implementation is similar to the implementation first disclosed in a NASA internal document [6] and U.S. Patent No. 6,895,115 B2. This implementation for two-dimensional data has been recently extended to accommodate one- and three-dimensional data [7]. This chapter is organized as follows. A high-level description of HSEG is followed by a more detailed description of RHSEG. Then a description of the parallel im-plementation of RHSEG is provided. Finally, timing comparisons are provided for several degrees of parallelism, from 256 CPUs down to 1 CPU. 5.2 Description of the Hierarchical Segmentation (HSEG) Algorithm The hierarchical image segmentation algorithm, HSEG, is based upon the relatively widelyutilizedhierarchicalstep-wiseoptimization(HSWO)regiongrowingapproach of Beaulieu and Goldberg [3], which can be summarized as follows: 1. Initialize the segmentation by assigning each image pixel a region label. If a pre-segmentation is provided, label each image pixel according to the pre-segmentation. Otherwise, label each image pixel as a separate region. 2. Calculatethedissimilaritycriterionvaluebetweenallpairsofspatiallyadjacent regions,findthepairofspatiallyadjacentregionswiththesmallestdissimilarity criterion value, and merge that pair of regions. 3. Stop if no more merges are required. Otherwise, return to step 2. HSEG differs from HSWO in one major aspect. The HSEG algorithm allows for the merging of spatially non-adjacent regions, as controlled by the Swght parameter. For Swght = 0.0, only spatially adjacent regions are allowed to merge, as in HSWO. However, for Swght > 0.0, HSEG allows merges between spatially non-adjacent regions. For Swght = 1.0, merges between spatially adjacent and non-adjacent re-gions are given equal weight. For values of Swght between 0.0 and 1.0, merges between spatially adjacent regions are favored by a factor of 1.0/Swght. Allowing for a range of merge priorities for spatially non-adjacent regions provides HSEG with a great deal of flexibility in tailoring the segmentation results to a particular need. © 2008 by Taylor & Francis Group, LLC 100 High-Performance Computing in Remote Sensing HSEG also provides a selection of dissimilarity functions for determining most similarpairsofregionsformerging.Theavailableselectionofdissimilarityfunctions is based on vector norms, mean-squared error, entropy, spectral information diver-gence (SID), spectral angle mapper (SAM), and normalized vector distance (NVD). SeetheRHSEGandHSEGViewerUser’sManual[8]forthemathematicaldefinitions ofthesedissimilarityfunctions.Optionsforotherdissimilarityfunctionscanbeeasily added. 5.3 The Recursive Formulation of HSEG Themergingofspatiallynon-adjacentregionsinHSEGleadstoheavycomputational demands. These demands are significantly reduced through a recursive approxima-tion of HSEG, called RHSEG, which recursively subdivides the imagery data into smaller sections to limit to a manageable number the number of regions considered at any point in the algorithm (usually in the range of 1000 to 4000 regions). RHSEG includes a provision to blend the results from the subsections to avoid processing window artifacts. This recursive approximation also leads to a very efficient paral-lel implementation. This parallel implementation of RHSEG is so efficient that full Landsat Thematic Mapper (TM) scenes (approximately 7000 by 6500 pixels) can be processed in 2 – 8 minutes on a Beowulf cluster consisting of 256 2.4GHz CPUs (http://thunderhead.gsfc.nasa.gov). This is only 10 to 20 times the amount of time that the Landsat TM sensor takes to collect this amount of data. The two spatial dimensional version of RHSEG was described in [5] and [9]. A description of RHSEG, generalized to ND spatial dimensions, follows (see also [7]): 1. Given an input image X, specify the number levels of recursion (Lr) required and pad the input image, if necessary, so that each spatial dimension of the data set can be evenly divided by 2(Lr −1). (A good value for Lr results in an image section at recursive level Lr consisting of roughly 1000 to 4000 pixels.) Set L = 1. 2. Call rhseg(L, X), where rhseg(L, X) is as follows: 2.1 If L = Lr, go to step 2.3. Otherwise, divide the image data into 2ND equal subsections and call rhseg(L + 1, X/2ND ) for each image section (represented as X/2ND ). 2.2 After all 2ND calls torhseg() from step 2.1 complete processing, reassem- ble the image segmentation results. 2.3 If L < Lr, initialize the segmentation with the reassembled segmentation resultsfromstep2.2.Otherwise,initializethesegmentationwithonepixel perregion.ExecutetheHSEGalgorithmontheimage X withthefollowing modification:Terminatethealgorithmwhenthenumberofregionsreaches the preset value Nmin. © 2008 by Taylor & Francis Group, LLC Parallel Implementation of the Recursive Approximation 101 3. Execute the HSEG algorithm (per Section 5.2) on the image X using as a pre-segmentation the segmentation output by the call to rhseg() in step 2. The defaults for the user specifiable parameters Lr and Nmin depend on the size of the image data and are calculated internally by RHSEG. Under a number of circumstances, the segmentations produced by the RHSEG algorithm exhibit processing window artifacts. These artifacts are region boundaries that are along the processing window seams, even though the image pixels across the seamsareverysimilar.Processingwindowartifactsareusuallyminorbutcanbemore noticeable, depending on the image. They tend to be more noticeable and prevalent in largerimages.However,allprocessingwindowartifactscanbecompletelyeliminated by adding a fourth step to the definition of rhseg(L, X) given above (following [10] and [11]): 2.4. If L = Lr, exit. Otherwise do the following (and then exit): (a) For each region, identify other regions that may contain pixels that are more similar to it than the region that they are currently in. These regions are placed in a candidate region label set for each region. This is done by: i. scanningtheprocessingwindowseambetweensectionsprocessed atthenextdeeperlevelofrecursionforpixelsthataremoresimilar (by a factor of Fseam) to the region existing across the processing window seam. ii. for Swght > 0.0 identifying regions that have a dissimilarity between each other less than F egion × Swght × Tmax (Tmax is the maximum of the merge threshold encountered so far in HSEG). (b) For each region with a candidate region label set of size greater than zero, identify pixels in the region that are more similar by a fac-tor of Fsplit to regions in the candidate region label set than to the region they are currently in. If Swght = 1.0, simply switch the region assignmentofthesepixelstothemoresimilarregion.Otherwise,split these pixels out of their current regions and remerge them through a restrictedversionofHSEGinwhichregiongrowingisperformedwith thesesplit-outpixelsandmergingisrestrictedtoneighboringregions, the region from which the pixel was split out from, and regions in the candidate region label set of the region from which the pixel was split out from. Processing window artifact elimination as introduced here not only eliminates the processing window artifacts, but does so with minimal computational overhead. Thecomputationtimenomorethandoublesforawiderangeofimagesizes[11].The default value of 1.5 for the parameters values Fseam, F egion, and Fsplit works well for a wide range of images. © 2008 by Taylor & Francis Group, LLC ... - tailieumienphi.vn