IMAGE segmentation and edge detection for multispectral

IEEE TRANSACTIONS ON IMAGE PROCESSING, VOL. 23, NO. 5, MAY 2014

2315

Model-Based Edge Detector for Spectral Imagery Using Sparse Spatiospectral Masks Biliana S. Paskaleva, Member, IEEE, Sebastián E. Godoy, Woo-Yong Jang, Steven C. Bender, Sanjay Krishna, Fellow, IEEE, and Majeed M. Hayat, Fellow, IEEE

Abstract— Two model-based algorithms for edge detection in spectral imagery are developed that specifically target capturing intrinsic features such as isoluminant edges that are characterized by a jump in color but not in intensity. Given prior knowledge of the classes of reflectance or emittance spectra associated with candidate objects in a scene, a small set of spectral-band ratios, which most profoundly identify the edge between each pair of materials, are selected to define a edge signature. The bands that form the edge signature are fed into a spatial mask, producing a sparse joint spatiospectral nonlinear operator. The first algorithm achieves edge detection for every material pair by matching the response of the operator at every pixel with the edge signature for the pair of materials. The second algorithm is a classifier-enhanced extension of the first algorithm that adaptively accentuates distinctive features before applying the spatiospectral operator. Both algorithms are extensively verified using spectral imagery from the airborne hyperspectral imager and from a dots-in-a-well midinfrared imager. In both cases, the multicolor gradient (MCG) and the hyperspectral/spatial detection of edges (HySPADE) edge detectors are used as a benchmark for comparison. The results demonstrate that the proposed algorithms outperform the MCG and HySPADE edge detectors in accuracy, especially when isoluminant edges are present. By requiring only a few bands as input to the spatiospectral operator, the algorithms enable significant levels of Manuscript received July 14, 2013; revised December 20, 2013 and February 27, 2014; accepted March 16, 2014. Date of publication April 1, 2014; date of current version April 22, 2014. This work was supported in part by the National Consortium for MASINT Research Partnership Project through the Los Alamos National Laboratory under Grant 57461-001-07, in part by the National Science Foundation under Grant ECS-0401154 and Grant IIS-0434102, in part by the Korea Research Institute of Standards and Science, in part by Nanostructures Far-IR/terahertz Detectors for Next Generation Imaging Technology Phase II Project, in part by the National Science Foundation under Grant ECCS-0925757, and in part by AFOSR Optoelectronic Research Center Grant. The work of S. E. Godoy was supported in part by the New Mexico Cancer Nanoscience and Microsystems Training Center, University of New Mexico, and in part by CONICYT, Chile. Sandia National Laboratories is a MultiProgram Laboratory managed and operated by Sandia Corporation, a wholly owned subsidiary of Lockheed Martin Corporation, through the U.S. Department of Energy National Nuclear Security Administration under Contract DE-AC04-94AL85000. The associate editor coordinating the review of this manuscript and approving it for publication was Dr. Debargha Mukherjee. B. S. Paskaleva was with the Center for High Technology Materials, Department of Electrical and Computer Engineering, University of New Mexico, Albuquerque, NM 87131 USA. She is currently with Sandia National Laboratories, Albuquerque, NM 87123 USA (e-mail: [email protected]). S. E. Godoy, S. C. Bender, S. Krishna, and M. M. Hayat are with the Center for High Technology Materials, Department of Electrical and Computer Engineering, University of New Mexico, Albuquerque, NM 87131 USA (e-mail: [email protected]; [email protected]; [email protected]; [email protected]). W-Y. Jang is with U.S. Air Force Research Laboratory, Wright-Patterson AFB, OH 45433 USA (email: [email protected]) Color versions of one or more of the figures in this paper are available online at http://ieeexplore.ieee.org. Digital Object Identifier 10.1109/TIP.2014.2315154

data compression in band selection. In the presented examples, the required operations per pixel are reduced by a factor of 71 with respect to those required by the MCG edge detector. Index Terms— Edge detection, isoluminant edge, classification, multicolor edge detection, spatio-spectral mask, spectral ratios.

I. I NTRODUCTION

I

MAGE segmentation and edge detection for multispectral (MS) and hyperspectral (HS) images can be an inherently difficult problem since gray-scale images associated with individual spectral bands may reveal different edges. Segmentation algorithms for gray-scale images utilize basic properties of intensity values such as discontinuity and similarity [1]. Popular gray-scale edge detectors include Canny [2], Sobel [3], and Prewitt [1], to name just a few. The transition from a gray-scale to a multicolor image complicates edge detection significantly: the standard definition of a gray-scale edge as a “ramp” or “ridge” between two regions [1, p. 573] is no longer appropriate because a multicolor image has multiple image planes (channels) corresponding to different spectral bands. Moreover, depending on the composition of the scene, two distinct spectral (color) regions may exhibit the same intensity for one or more bands and, in this case, the edge between the two regions is termed isoluminant. An isoluminant edge is therefore characterized by a jump in color rather than a jump in intensity. As a result, isoluminant edges cannot be detected easily by a standard gradient-based operator because they usually do not exhibit an intensity ramp that can be estimated by the magnitude of such an operator [4]. (Examples of isoluminant edges will be shown in Section III-B). The extension of gray-scale edge detection to multicolor images has followed three principal paths [5]. A straightforward approach is to apply differential operators, such as the gradient, separately to each image plane and then consolidate the information to obtain edge information [6]. Sapiro [7], Sandberg and Chan [8] identified several key drawbacks of such a straightforward approach. First, edges can be defined by combinations of different image planes and these edges may be missing in some of the image planes. Second, processing image planes separately disregards potential correlation across image planes. Third, integration of information from separate image planes is not trivial and is often done in an ad hoc manner. Moreover, in cases when an edge appears only in a subset of image planes, there are no standard ways to fuse the information from different planes. In recent years, new gray-scale algorithms were presented in the community (see, for example, [9]–[12]) but they all

1057-7149 © 2014 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission. See http://www.ieee.org/publications_standards/publications/rights/index.html for more information.

2316

suffer from the scalability problem when applied to multicolor images. A second approach for multicolor edge detection is to embed the variations of all color channels in a single measure, which is then used to obtain the edge maps [5]. Typically, this approach is developed by starting from a given gray-scale operator, which is then consistently extended to multicolor images. Two representative examples of this approach are the multicolor gradient (MCG), proposed by Di Zenzo [13], and the morphological color gradient (MoCG) of Evans and Liu [14]. The MCG operator represents a consistent extension of the standard gradient operator to multicolor images and it measures the local steepness of the multicolor image considered as a manifold embedded in the Euclidean space. (A hyper-pixel belongs to a multicolor edge if the local steepness of the manifold, as measured by MCG, exceeds a given threshold [13]). Similarly, the MoCG operator is a consistent extension of the morphological gray-scale gradient operator to multicolor images [14], [15]. Such an operator is defined as the difference of the dilation and the erosion operators applied to a given structuring element [16]. Because the MCG and the MoCG edge detectors simultaneously utilize spatial and spectral information, they are examples of joint spatio-spectral image-processing algorithms. The MCG algorithm and its related algorithms have been used with great success in digital image processing [7], [17], [18]; however, as shown by the complexity estimates in Section IV, for spectral images with a large number of bands, the number of operations required by the MCG algorithm can be prohibitively high. Another approach that falls into the category of joint spatiospectral algorithms is the order-statistics approach [19] and its extensions presented by Toivanen et al. [20] and Jordan and Angelopoulou [21]. In general, these algorithms consider the data as a discrete vector field and they utilize an R-ordering method to define a color edge detector using the magnitudes of linear combinations of the sorted vectors. Another algorithm, presented by Yu et al. [22], performs a mapping of the color image into a feature space, which considers features such as the local contrast, the edge connectivity, the color contrast similarity and the orientation consistence. These features are merged together to create a single feature that is compared to a threshold to generate the final edge map. Similar approaches can be found in [23] and [24] and more sophisticated estimators of the color gradient can be found in [25]. A third approach for multicolor edge detection is to aim the algorithm to detect solely the changes between the materials present in the imaged scene. For example, the HySPADE (Hyperspectral/Spatial Detection of Edges) algorithm transforms the data cube into a spectral angle (SA) cube by calculating the SA between each hyperpixel in the cube with every other pixel [26]. As a result, the third dimension in the original cube is replaced by the SA results, where jumps in the SA represent changes in the materials. The positions of these jumps are mapped back in to the original data cube, and the final edge map is derived upon statistical accumulation of edge information contained in every SA-cube. One important distinction between the MCG algorithm and the HySPADE algorithm is that the former utilizes both spectral and spatial

IEEE TRANSACTIONS ON IMAGE PROCESSING, VOL. 23, NO. 5, MAY 2014

information to detect the edges, while the HySPADE utilizes solely spectral information to unveil the boundaries of the material composition. In this paper, we propose two novel joint spatio-spectral algorithms based on the concept of spectral-bands ratios. By judiciously selecting the few and most relevant spectral bands to maximize the contrast between pairs of materials, we define a small set of ratios of spectral-band features that most profoundly identify edges between each pair of materials. In order to obtain the collection of relevant bands and the corresponding ratios, a collection that we term the edge signature, we utilize prior knowledge of the spectral characteristics of the materials in the scene, obtained, for example, from a library of spectral data. Through this stage, the proposed algorithms achieve substantial levels of data compression as compared with the MCG or HySPADE algorithms. In conjunction with a spatial mask, the few spectral bands from the edge signature give rise to a multispectral operator that can be viewed as a sparse, three-dimensional (3D) mask, which is at the heart of the two proposed edge-detection algorithms. The 3D mask does not operate on a single image plane but instead it fuses spectral information from multiple image planes with spatial information. The first algorithm, termed the Spectral Ratio Contrast (SRC) edge detection algorithm, defines the edge map of a spectral image by matching the output of the 3D mask with the ratios from the edge signature. The second algorithm is an extension of the SRC algorithm and utilizes spectral classification to further enhance the detection of edges that are solely due to material (not intensity) changes. We term this extension the Adaptive Spectral Ratio Contrast (ASRC) edge detection algorithm since it adaptively changes the SRC algorithm sensitivity to edges (at each pixel) by considering the material-classification results of the neighboring pixels. The first algorithm detects edges that arise from both intensity and spectral changes while the second algorithm detects edges based on spectral changes only. The sparse spatio-spectral mask used in the SRC and the ASRC algorithms is an important mark of distinction from the MCG-based edge detector and other multispectral edgedetection algorithms. A second key distinctive mark of the proposed two algorithms is that they are not derivative-based: edge detection is effected by matching an edge signature rather than by measuring the gradient’s magnitude. Moreover, the application of spectral ratios to define multispectral operators for edge detection is a novel and a previously unexplored research direction. However, spectral ratios have been previously used in many techniques for quantitative vegetation monitoring [27]–[29], regional seismic discrimination [30]–[33] and deblurring of noisy multichannel images [34]. The paper is organized as follows. In Section II we present the SRC and the ASRC algorithms. In Section III we present results of applying the algorithm to real data from the Airborne Hyperspectral Imager (AHI) [35] and a quantum dots-in-awell (DWELL) mid-infrared (IR) imager [36] and compare the performance to those resulting from the Canny, MCG and HySPADE edge detectors. In Section IV we present a complexity analysis of the algorithms. Our conclusions are presented in Section V.

PASKALEVA et al.: MODEL-BASED EDGE DETECTOR FOR SPECTRAL IMAGERY

II. S PECTRAL R ATIO C ONTRAST A LGORITHM FOR E DGE D ETECTION An MS or HS image, also termed an image cube, is a 3D array of real numbers that we denote by u ∈ IR I ×J ×K , where I and J represent the number of horizontal and vertical pixels, respectively, in the spatial domain, and K represents the number of spectral bands. We denote any element of the cube u as u k (i, j ), where 1 ≤ i ≤ I , 1 ≤ j ≤ J and 1 < k ≤ K . The value of u k (i, j ) is referred to as the intensity of the kth band at spatial location (i, j ). For a fixed spatial location (i, j ), the K -dimensional vector u(i, j ) = (u 1 (i, j ), . . . , u K (i, j )) is termed a hyper-pixel. Meanwhile, for a fixed band index k, the two-dimensional array u k (·, ·) defines the kth image plane (color slice) of the spectral image. The goal of this paper is to define an edge map, F : IR I ×J ×K → {0, 1} I ×J, that assigns the value 1 to the pixel location (i, j ) if u(i, j ) belongs to an edge, while assigning the value 0 otherwise. In this paper we define an edge to be either a jump in the broadband intensity (as in the conventional definition of an edge for gray-scale images) or a change in the material that exhibit color contrast but not necessarily luminance contrast (as in isoluminant edges). The SRC algorithm is designed with the objective of detecting both types of edges. To reduce the detection of false edges triggered by noise, and to improve the detection of isoluminant edges, the ASRC algorithm is designed to detect edges that are purely due to changes in the material. For two types of materials A and B, we build the edge map FAB in three stages: (i) model-based edge signature identification, (ii) sparse spatio-spectral mask development, and (iii) edge discrimination. The SRC and the ASRC algorithms share the first two stages, but they differ in the third stage, where the latter involves spectral classification before discriminate the edges. For simplicity, we first describe both algorithms assuming only two distinct materials A and B in scenes, and later we describe the extensions to multiple materials. A. Model-Based Edge Signature Identification Given two distinct materials A and B in a scene that is probed by a sensor, we seek those bands for A and B whose ratios can best identify the spatial locations that correspond to the boundary points between the two materials. Let E AB denote the set of all spatial indices (in the cube resulting from imaging the scene containing materials A and B) that are boundary points between materials A and B, and let a = (a1 , . . . , a K ) and b = (b1 , . . . , b K ) denote hyper-pixels corresponding to materials A and B, respectively. For example, the vector a can be obtained by taking the average of all hyperpixels as material A is probed by the sensor. We next define the spectral ratio index between materials A and B as the K × K matrix ⎛ a1 a1 ⎞ ··· bK ⎟ ⎜ b1 ⎜ .. .. ⎟. . . (1) A/B  ⎜ . . . ⎟ ⎝a aK ⎠ K ··· b1 bK

2317

We define the signature of the edge between materials A and B as a small collection of size R, where R  K , of elements of (1) that can reliably identify, as described below, the edge’s spatial indices, E AB . We denote the edge signature by EAB = {( p1, q1 , ρ1 ), . . . , ( p R , q R , ρ R )},

(2)

where pr and qr are the band indices associated with the ratios ρr = a pr /bqr , r = 1, . . . , R. The integer R is the length of the edge signature. The selection of the triplets ( pr , qr , ρr ), r = 1, . . . , R, is made as follows. First, we select S bands {i 1 , . . . , i S }, where the materials A and B exhibit maximum separation, i.e., i 1 = arg max |ai − bi |, 1≤i≤K

i 2 = arg max |ai − bi |, 1≤i≤K , i =i1

.. . i S = arg

max

1≤i≤K i =i1 ,...,i S−1

|ai − bi |.

Clearly, in the case of HS imagery the search for the best S bands will be in a larger space compared to the case of MS imagery, but the same procedure is utilize in both types of imagery. Hence, once the edge signature is obtained the complexity associated with processing MS and HS imagery are identical. Next, we compute the spectral ratios using all possible band combinations: ρ pq = ai p /biq , 1 ≤ p, q ≤ S . Without loss of generality, we may assume that all ratios are less than or equal to unity. (If ρ pq > 1 for some p and q, we simply replace it by its reciprocal.) Finally, we define EAB by selecting the R ratios that exhibit the strongest spectral contrast between the classes. To rank the ratios according to their spectral contrast, we note that owing to the convention that ρ pq ≤ 1, the ratios closest to zero correspond to the strongest spectral contrast between any two bands. Thus, we select the first pair of bands, { p1, q1 }, as the pair corresponding to the smallest ratio, ρ1 = ρ p1 q1 = arg min ρ pq , 1≤ p,q≤S

the second pair of bands { p2 , q2 } as the pair corresponding to the next smallest ratio, ρ2 = ρ p2 q2 = arg min ρ pq , 1≤ p,q≤S p = p1 q =q1

and so on. We combine the ordered band indices and the corresponding ranked ratios to define the R triplets in the edge signature (2). Since all the ratios are less than or equal to unity, it can be easily shown that the definition of the edge signature is invariant under the change in the order of the materials A and B. B. Sparse Spatio-Spectral Mask Development We denote a spatial mask, M, at a pixel (i, j ) as a list of pixel pairs surrounding the pixel of interest. More precisely, M(i, j ) is the union of M neighborhoods of

2318

IEEE TRANSACTIONS ON IMAGE PROCESSING, VOL. 23, NO. 5, MAY 2014

M pixels, M(i, j ) = ∪m=1 Nm (i, j ), where each neighborhood consists of two distinct pixels surrounding (i, j ), i.e., Nm (i, j ) = {u m− (i, j ), u m+ (i, j )}. For example, one can define a 3 × 3 mask centered at the pixel (i, j ) that excludes the center pixel by taking the union of four neighborhoods, N1 = {u 1− (i, j ), u 1+ (i, j )} = {u(i − 1, j ), u(i + 1, j )}, N2 = {u 2− (i, j ), u 2+ (i, j )} = {u(i, j − 1), u(i, j + 1)}, etc. Next, we define the operation of the joint spatio-spectral mask at the (i, j )-th pixel by computing the ratios between each of the M pixel pairs of the spatial mask M, at each of the R band pairs given in the edge signature E AB . For example, using the pair of hyper-pixels defined by N1 and the pair of bands given by first triplet of the edge signature, 1+ ( p1 , q1 , ρ1 ), one can define the ratios u 1− p1 (i, j )/u q1 (i, j ) and 1+ 1− u p1 (i, j )/u q1 (i, j ). Namely, the application of the spatiospectral mask to each location (i, j ) results in a 2M × R matrix of “features.” (Compare this to a gray-scale image when the application of a spatial mask to a pixel results in a scalar). Now the application of the spatio-spectral mask to the entire image cube defines the mapping  I ×J , (3) KAB : IR I ×J ×K → IR2M×R

where the (i, j )th entry of KAB (u) matrix of spectral ratios given by ⎛ u 1− p1 (i, j ) ⎜ 1+ ⎜ u q1 (i, j ) ⎜ ⎜ .. ⎜ ⎜ M−. ⎜ u (i, j ) ⎜ p1 ⎜ M+ ⎜ u q1 (i, j ) KAB (u)(i, j ) = ⎜ ⎜ u 1+ (i, j ) ⎜ p1 ⎜ 1− ⎜ u q (i, j ) ⎜ 1 ⎜ .. ⎜ . ⎜ ⎜ u M+ (i, j ) ⎝ p1 u qM− 1 (i, j )

will be a 2M × R feature

··· ..

.

··· ··· ..

.

···

u 1− p R (i, j )

⎞

⎟ u 1+ q R (i, j ) ⎟ ⎟ .. ⎟ ⎟ . ⎟ ⎟ (i, j ) u M− pR ⎟ ⎟ M+ u q R (i, j ) ⎟ ⎟. ⎟ u 1+ p R (i, j ) ⎟ ⎟ ⎟ u 1− q R (i, j ) ⎟ ⎟ .. ⎟ . ⎟ M+ u p R (i, j ) ⎟ ⎠ u qM− (i, j ) R

Fig. 1. Representative joint spatio-spectral mask, KAB (u)(i, j), for the neighborhoods N1 and N2 described in the text and edge signature EAB of length 1 (M = 2 and R = 1). We utilize bold arrows to represent the first two ratios and regular arrows to represent the last two ratios in the presented example.

To illustrate the functionality of the joint spatio-spectral mapping KAB , consider a simple example for which the spatial mask is comprised of the simple neighborhoods N1 and N2 described earlier and an edge signature whose length is unity (M = 2 and R = 1). Note that in this example KAB (u)(i, j ) is a 4 × 1 matrix given by ⎛ 1− ⎞ u p1 (i, j )/u 1+ q1 (i, j ) ⎜ ⎟ ⎜ u 2− (i, j )/u 2+ (i, j ) ⎟ ⎜ p1 ⎟ q1 ⎟. (5) KAB (u)(i, j ) = ⎜ ⎜ 1+ ⎟ ⎜ u p1 (i, j )/u 1− ⎟ (i, j ) q 1 ⎝ ⎠ 2− u 2+ p1 (i, j )/u q1 (i, j )

(4)

For convenience, we denote the entries of the matrix KAB (u)(i, j ) as m+ κAB (i, j ; m − , r ) = u m− pr (i, j )/u qr (i, j ),

A visual illustration of the computation of the rows elements in this matrix is shown in Fig. 1. This figure shows that, unlike a conventional mask, KAB (u) does not operate on a single image plane of the image cube but instead it fuses information from different planes in a nonlinear fashion. Note that the spatiospectral mask is generally non-separable, i.e., it cannot be written as the product of an operator acting solely in the spatial domain and another operator acting in the spectral domain. Next, we develop the edge discrimination stage for the SRC algorithm. Later, the same stage is modified to generate the ASRC algorithm.

and m− κAB (i, j ; m + , r ) = u m+ pr (i, j )/u qr (i, j ).

The use of both κAB (i, j ; m − , r ) and κAB (i, j ; m + , r ) in (4) is required to account for the two possible material configurations at the mth hyper-pixel pair. Specifically, the first ratio captures the case when the hyper-pixel u m− (i, j ) is, for example, from material A and u m+ (i, j ) is from material B, whereas the second ratio is needed to account for the possibility that u m− (i, j ) is from material B and u m+ (i, j ) is from material A. Therefore, the use of the two ratios removes dependence on the direction of the transition between A and B, and it is similar to the use of the magnitude in the gradient operator to achieve its rotational invariance.

C. Edge Identification The third stage of the SRC algorithm is the utilization of the sparse spatio-spectral mask, KAB , to identify the edges between materials A and B. The proposed edge-identification process is based on the following rationale. In the ideal case when no noise is present and the image under test is comprised only of hyper-pixels with the exact same value of the characteristic hyper-pixels a and b, the output of the spatio-spectral mask will perfectly match the values of the ratios obtained from the edge signature. To illustrate this point, assume the same example given earlier, where the spatiospectral mask is given by (5). When we have a horizontal edge, the second and fourth ratios will have a value that is not

PASKALEVA et al.: MODEL-BASED EDGE DETECTOR FOR SPECTRAL IMAGERY

meaningful, but the first or third ratio will match the ratio from the edge signature, ρ p1 q1 = a p1 /bq1 . Indeed, if the upper pixel is from material A and the lower pixel is from material B, 1+ then the first entry of (5) will be u 1− p1 (i, j )/u q1 (i, j ) = a p1 /bq1 (which matches ρ p1 q1 ) and the third entry will be 1− u 1+ p1 (i, j )/u q1 (i, j ) = b p1 /aq1 . Conversely, when the upper pixel is from material B and the lower pixel is from material A, the computed ratios are switched, which means that the third entry of (5) will match ρ p1 q1 . This example shows that if an edge is present then at least one row of KAB (u)(i, j ) will perfectly match the ratios from the edge signature. When noise is present, we allow a matching tolerance to account for the similarity between the outcome of the mask and the edge signature ratios. This can be accomplished by defining the mapping I ×J I ×J   → {0, 1} M×R ,  : IR2M×R where the (i, j )th entry of KAB (u) is binary indicator matrix ⎛ 1 δ1 (i, j ) ⎜ δ 2 (i, j ) ⎜ 1 

⎜  KAB (u) (i, j ) = ⎜ .. ⎜ . ⎝ δ1M (i, j )

and the entries δrm (i, according to the rule ⎧ 1, ⎪ ⎨ m δr (i, j ) = 1, ⎪ ⎩ 0,

used to form the M × R ··· ··· ..

.

···

δ 1R (i, j )

⎞

δ 2R (i, j ) ⎟ ⎟ ⎟ ⎟ .. ⎟ . ⎠ δ RM (i, j )

(6)

j ) are assigned the values of 0 or 1 if |κAB (i, j ; m − , r ) − ρr | < , if |κAB (i, j ; m + , r ) − ρr | < ,

(7)

otherwise.

Here, the tolerance parameter, , accounts for both the natural variability and the presence of noise in the spectral data for materials A and B. Ideally, if the mth hyper-pixel pair belongs to the same material type, then the test in (7) will return the value of zero. Conversely, if the hyper-pixels forming the pair are from the two different materials, either the entry κAB (i, j ; m − , r ) or κAB (i, j ; m +, r ) will be equal to the corresponding ratio ρr from the edge signature EAB . As a result, the above test will return the value 1 for the elements δrm in the mth row of (6). As such, for a given pair of pixels, the number of non-zero elements in the associated mth row of the indicator matrix reveals the number of times the response of the mask KAB (u)(i, j ) has matched (within the specified tolerance ) the spectral ratios from the edge signature EAB . Because the pixel pairs used to form the rows of the mask correspond to different edge orientations (horizontal, vertical or diagonal), the number of ones in each row of (6) indicates the strength of the edge at position (i, j ) for that particular direction. One way to account for such strength is by computing the matrix infinity norm of (6). Specifically, we define the mapping I ×J  → {0, 1} I ×J ,  : {0, 1} M×R

2319

which converts the indicator matrix (6) into an edge map by 

  

1, if  KAB (u) (i, j ) ∞ ≥ R    KAB (u) (i, j ) = 0, otherwise, (8) R where for any matrix A, A ∞ = max1≤i≤M j =1 |ai j |, and  ≤ R is a specified threshold. With the definition in (8), the R (i, j ) location will belong to the collection E AB of edges if the edge strength in at least one direction, as measured by the  number of ones in the rows of (6), exceeds the threshold R.  then the (i, j ) site If none of the edge strengths exceeds R,  can be does not belong to E AB . The value of the threshold R used to adjust the sensitivity of the edge detector to noise. For  makes the algorithm less sensitive to example, increasing R noise but more restrictive. Finally, we define the edge map as the composition FAB =  ◦  ◦ KAB . Note that KAB is the only problem-specific component in FAB ; the functions  and  are not problem specific. Because the edge signatures are determined independently for each pair of materials and the information from different color slices is properly fused, the SRC algorithm is particularly well-suited for scenes that contain both weak edges (e.g., isoluminant edges) and strong edges. In contrast, an MCGbased edge detector would require, for the detection of weak edges, reducing the threshold at the expense of producing false edges, as discussed in Section III. In the next subsection we present an extension of the edge identification stage of the SRC algorithm to benefit from locally-adaptive thresholding based on material classification. D. Classifier-Enhanced Edge Discrimination As shown in Section III, the SRC algorithm is capable of detecting edges that are due to intensity changes and edges that are due to spectral changes. In this subsection, we develop the ASRC algorithm, which is restricted to capturing edges that are due solely to material changes (and not intensity changes). By specializing the algorithm to material changes, we increase its tolerance to noise and the corresponding false edges are minimized. In order to capture solely the changes between materials in the ASRC algorithm, we utilize material classification of neighboring pixels to adaptively pre-qualify the spectral ratios before computing the indicator matrix (6). The use of classification to enhance other tasks such as segmentation is an area that has already been studied [37], [38]. For example, Loog and Ginneken [37] utilized a k-nearest neighbor classifier to generate an initial segmentation of ribs in chest radiographs, which is iteratively updated using other features such as spatial distribution of the pixels by means of different classifiers. Here, we fuse classification into the SRC algorithm to minimize the effect of the misidentified pixels from the former algorithm and improve the edge identification stage of the latter. Let us consider the feature matrix of the (i, j )th location, KAB (u)(i, j ), as given in (4). We want to rank the entries

2320

IEEE TRANSACTIONS ON IMAGE PROCESSING, VOL. 23, NO. 5, MAY 2014

κAB (i, j ; m ∓, r ) of the feature matrix with the following two objectives in mind: (i) promoting the thresholding of ratios (as edge candidates) when the spatio-spectral mask contains hyper-pixels from two distinct materials and (ii) discouraging the thresholding of ratios when the mask contains hyper-pixels from only one type of material. To do so, we embed the i, j data-dependent multiplicative factors γm,r with the tolerance  in (7), which results in its redefinition as ⎧ i, j ⎪ 1, if |κAB (i, j ; m − , r ) − ρr | < γm,r , ⎪ ⎨ i, j δrm (i, j ) = 1, if |κAB (i, j ; m + , r ) − ρr | < γm,r , ⎪ ⎪ ⎩ 0, otherwise, i, j

where we use the same multiplicative factor, γm,r , for both tests in order to maintain the independence on the direction of the transition between materials. The mathematical i, j definition of the multiplicative factors γm,r is as follows. Let Nu (i, j ), Nl (i, j ), Nr (i, j ) and N (i, j ) denote four neighborhoods surrounding the (i, j )th hyper-pixel of interest. The subscripts u, l, r and stand for upper, lower, right and left neighborhoods, respectively. The understanding is that Nu (i, j ) contains neighboring pixels above the hyperpixel (i, j ), N (i, j ) contains neighboring pixels to the left of pixel (i, j ), and so on. We do not impose any restrictions on the four neighborhoods at this point. Let the function C : IR S → {0, 1} be a classifier that maps each hyper-pixel to a class of materials, where outputs 0 and 1 represent classes A and B, respectively. Now for any collection of indices N , we denote the class of N by C(N ), which is defined according to a certain prescribed classification rule. Next, define   i, j C(Nr ) ⊕ C(N ) , (9) γm,r = C(Nu ) ⊕ C(Nl ) where  ⊕ denotes the “exclusive OR” operation and the symbol “ ” represents the “OR” operation. For simplicity of the notation, we have discarded the (i, j ) dependence of each neighborhood set N in (9) with the understanding that each neighborhood is defined on a pixel-by-pixel basis. i, j From (9), γm,r will be unity (in which case the (m, r )th pixel-band pair at the (i, j ) location qualifies for thresholding as usual) if at least one of the opposite neighborhoods are i, j classified as two different materials. On the other hand, γm,r will be zero (in which case the (m, r )th pixel-band pair at the (i, j ) location does not qualify for thresholding) if the declared class of each neighborhood is in agreement with the declared class of its opposite neighborhood. As a consequence, the ASRC will operate as the the SRC algorithm only if the outcome of the classifier indicates the possible presence of an edge, suppressing edges that are due to an intensity change. This will also reduce the detection of false edges. The fact that the ASRC algorithm is restricted to identifying edges based on color only is similar to that of the HySPADE algorithm [26]; however, the algorithms are conceptually different. A key difference between the ASRC and the HySPADE algorithms is that the former utilizes the sparse, 3D mask of ratios to fuse spectral and spatial information to nonlinearly extract edge information while the latter algorithm utilizes only spectral information to compute spectral

angles, which are linear spatio-spectral features. Even though the HySPADE algorithm performs equally well compared to the SRC algorithm, its performance is worse than that of the ASRC algorithm in the presence of isoluminant edges (see, for example, Fig. 6). Moreover, as presented in the Section IV, the HySPADE algorithm requires a high number of operations per pixel (>109 operations per pixel in our examples), as compared with the operations required by the proposed algorithms (< 90 operations per pixel in our examples). E. Extension of the Algorithms to Multiple Materials We limit the description of the extension of the algorithm to multiple materials to the case of three distinct materials A, B and C. The extension to the general case is straightforward. Due to the invariability of the detector for the order of the materials, for three distinct materials A, B and C there are three possible edges: E AB , E AC and E BC . In this case we obtain three edge signatures, EAB , E AC and E BC , from which we define three joint spatio-spectral masks KAB , K AC , and K BC . We then use these masks to identify the hyper-pixels belonging to the edge E AB between materials A and B, the hyper-pixels from the edge E AC between materials A and C, and the hyper-pixels from the edge E BC between materials B and C. The final edge map is obtained by the union of the three edges: E ABC = E AB ∪ E AC ∪ E BC .

III. E XPERIMENTAL R ESULTS In our study, we employ raw HS imagery from the AHI sensor, and raw MS imagery from the DWELL sensor (more details to follow). In order to create a more challenging scenario for the algorithms, we normalize the data by their broadband intensity [39]. The normalization minimizes the role of broadband emissivity in the discrimination process and emphasizes the spectral contrast. For the AHI dataset, we only perform a qualitative comparison of the algorithms since the ground-truth information is not available for this dataset. On the other hand, for the data from the DWELL sensor we perform both qualitative and quantitative assessment of the proposed algorithms and the benchmark algorithms as the ground-truth information is available. We compare the outcome of our algorithms with the edge maps obtained by the Canny algorithm [2] (applied to selected bands), the MCG algorithm [13] and the HySPADE algorithm [26]. We restrict our attention to edge signatures with unity length using two bands (i.e., S = 2 and R = 1), which is the minimum value required by the algorithms. Moreover, we utilize a 3 × 3 spatial mask to construct the joint spatio-spectral mask, KAB . Within the spatial mask, we identify four directions (each one associated with a pair of pixels): horizontal, vertical and the two diagonals, i.e., M = 4. For the ASRC algorithm, we select the distance-based Euclidean classifier for its simplicity and the good results observed; the neighborhood sets Nu (i, j ), Nl (i, j ), Nr (i, j )

PASKALEVA et al.: MODEL-BASED EDGE DETECTOR FOR SPECTRAL IMAGERY

2321

TABLE I T HE E DGE S IGNATURES B ETWEEN C LASSES B, G AND R O BTAINED FOR THE AHI D ATA

and N (i, j ) are defined within the same 3 × 3 spatial mask used in the SRC algorithm. This choice of spatial mask, classifier and neighborhood sets is also considered for the complexity analysis in Section IV. A. Edge Detection Using AHI Imagery The AHI sensor consists of a long-wave IR (7μm-11.5μm) pushbroom HS imager and a visible high-resolution CCD linescan camera. The HS imager has a focal-plane array (FPA) of 256 × 256 elements with spectral resolution of 0.1μm [35]. For this study, we utilize AHI data that contains three different classes: building (B), ground (G) and road (R). We utilize the 200 low-noise bands out of the 256 available bands. The calculated edge signatures (band indices and the corresponding responding ratios) for each pair of materials are summarized in Table I. Figure 2 shows a comparison among the edge maps obtained by the Canny, the HySPADE, the MCG, the SRC and the ASRC algorithms for the raw sensor data (first and second columns) and for the normalized data (third and fourth columns). The Canny algorithm is applied to the same depicted image, which corresponds to the image plane at band 14. Recall that the MCG and the SRC algorithms detect edges characterized by both intensity and spectral changes. The HySPADE and ASRC algorithms, on the other hand, detect edges that exhibit a change in the spectral content only. From the results presented in Fig. 2 we observe that the Canny edge detector performs very well when applied to the AHI raw image for spectral band 14 (row I, column b). However, when the algorithm is applied to the intensitynormalized image, the performance of the Canny algorithm significantly degrades (row I, column d). This degradation is a result of the fact that the Canny algorithm detects intensity changes only, and it is expected to perform optimally for high intensity contrast images such as the image in row I, column a. The MCG and the SRC algorithms produce virtually the same edge maps when applied to raw sensor data (second row, columns a and b), with a clear computational advantage seen in the SRC algorithm by requiring only two spectral bands, whereas the MCG algorithm requires all the 200 available bands. When normalized data is used (second row, c and d columns), few edges in some areas are missed either by the SRC or the MCG algorithms. Nonetheless, the edge maps between the two algorithms are again comparable. Moreover, the results for the normalized case are very similar to those

Fig. 2. Comparison between the SRC and ASRC edge detectors with the edge signatures defined in Table I on the one hand, and the Canny, MCG and HySPADE edge detectors on the other hand, for the raw sensor data (columns a and b) and the normalized sensor data (columns c and d) for the AHI imagery. First row, columns a to d: raw AHI data (band 14); Canny edge map for raw AHI data at band 14; normalized AHI data at band 14; Canny edge map for AHI normalized data at band 14; Second row, columns a to d: MCG and SRC edge maps for AHI raw dataset; MCG and SRC edge maps for AHI normalized dataset; Third row, columns a to d: HySPADE and ASRC edge maps for AHI raw dataset; HySPADE and ASRC edge maps for AHI normalized dataset.

for the raw data case. These results show the advantage of the methods that utilize both intensity and spectral information over purely gray-scale algorithms such as Canny.

2322

The ASRC algorithm (row III, column b) performs significantly better compared to the HySPADE algorithm (row III, column a) when applied to the AHI raw data. The edge map obtained by HySPADE exhibits noise and some of the edges that were detected by SRC, MCG and ASRC are missed by HySPADE. The advantage of ASRC over HySPADE continues to be pronounced when both algorithms are applied to the normalized data (row III, columns d and c). The edge maps obtained by the ASRC algorithm applied to raw and normalized data (row III, columns b and d) are virtually identical. This is due to the fact that the ASRC algorithm detects edges based on changes of the spectral content only. As for HySPADE, the application of the algorithm to the normalized AHI data results in a slight degradation of the edge detection but overall reduction of the noise in the edge map compared to application to the raw AHI data (row III, column a); however, as in the case of ASRC, the edge maps are comparable. It is important to note that the edge maps obtained by the MCG, the SRC and the ASRC algorithms are very similar for both raw and normalized AHI data cases. One important conclusion can be drawn from the results presented so far. By choosing only a few bands with maximum spectral separation and by allowing unrestricted band combinations to form the ratios, the SRC and ASRC algorithms (with edge signatures that use the minimal possible length) perform as well as the MCG algorithm and outperform the HySPADE algorithm. This is an important result because it lends itself to substantial data compression, compared to MCG, as well as fast processing, as compared to HySPADE. The SRC and ASRC algorithms offer a performance advantage over the Canny, the MCG and the HySPADE algorithms for images that contain isoluminant edges as seen next for the DWELL imagery. B. Edge Detection Using DWELL Imagery The DWELL sensor used in these experiments was designed and fabricated at the Center for High Technology Materials at the University of New Mexico [36], [39]. The DWELL photodetector offers a unique property of spectral tunability that is continuously controllable through the applied bias voltage. This feature of the DWELL is a result of the quantum-confined Stark effect [40]. In essence, a single DWELL photodetector can be thought of as a continuously tunable MS spectral detector, albeit with overlapping spectral bands [39]. In these experiments we utilize a 320×256 DWELL FPA to image two different arrangements of rocks, as shown in Fig. 3 (first column). The first arrangement (top-left) is comprised of granite (G) and limestone (L) rocks (approximately 1–2 inch in diameter). The surrounding background (B) in this image corresponds to the opening of a blackbody source. The second arrangement (bottom-left) is comprised of the rocks phyllite (P), granite (G) and limestone (L), surrounded by the same background (B) as that in the first arrangement. Both examples contain an invisible isoluminant edge between the granite and the limestone rocks that exists on the tip of the black arrows. The edge maps shown in Fig. 3 were obtained by using the Sobel (second column) and the Canny (third column)

IEEE TRANSACTIONS ON IMAGE PROCESSING, VOL. 23, NO. 5, MAY 2014

Fig. 3. Two datasets used in the current study: the first dataset is comprised of B, G and L classes (top row) and the second dataset is comprised of B, G, L and P classes (bottom row). First column: images acquired with the DWELL FPA (with enhanced contrast to show details) operating at an applied bias of 1.0V . The isoluminant edges (not visible) are marked by the tips of the black arrows. Second column: edge map obtained by the Sobel gray-scale edge detector; third column: edge map obtained by the Canny gray-scale edge detector.

Fig. 4. Spectral response of the DWELL photodetector at an applied bias of 1.0 V.

edge detectors applied to raw DWELL-sensor data when the FPA is operated at 1.0V . The corresponding spectral response of the sensor at the applied bias of 1.0V is shown in Fig. 4. Note that the Sobel edge detector has entirely missed the edge between granite and limestone rocks in both examples. Moreover, it has also failed to detect strong edges between both the granite-phyllite pair and the limestone-phyllite pair. However, the more sophisticated Canny edge detector picks up these strong edges, and it partially detects the isoluminant edge in the first examples. Nevertheless, it does not detect the isoluminant edge in the second example. By operating the DWELL sensor at ten different bias voltages, we generated a multispectral cube to test the proposed algorithms. The obtained edge signature triplets for all the possible combinations of material pairs for both datasets are summarized in Table II. In what follows, we will term the DWELL imagery that contains background, granite and limestone classes, as shown in Fig. 3 (top-left), the first DWELL dataset, and we term the imagery that contains background, phyllite, granite and limestone, as shown in Fig. 3 (bottomleft), the second DWELL dataset. The results for the first DWELL dataset for raw sensor data are shown in Fig. 5. The first row of edge maps shows the results of the application of the Canny edge detector to four randomly selected bands. It is important to note that some bands present a high number of false edges, whereas for other bands the isoluminant edges are detected. As such,

PASKALEVA et al.: MODEL-BASED EDGE DETECTOR FOR SPECTRAL IMAGERY

2323

TABLE II E DGE S IGNATURES A MONG THE B, P, G, AND L C LASSES

TABLE III C OMPARISON TABLE FOR THE PD AND PF R ESULTS OF F IVE

O BTAINED FOR THE DWELL D ATASETS

A LGORITHMS (C ANNY, MCG, H Y SPADE, SRC AND ASRC) FOR THE

Fig. 5. Comparison between the Canny algorithm applied to individual bands (first row), MCG algorithm (second row), HySPADE algorithm (third row), SRC algorithm (fourth row) and ASRC algorithm (fifth row) for the dataset containing granite and limestone rocks (first dataset). The Canny algorithm was applied to the images at bands 1, 6, 8 and 9, respectively. The MCG and HySPADE results are presented for a sequence of increasingly permissive tolerances in order to unveil the isoluminant edge. Last two rows show the SRC and ASRC edge maps: first column, the edges E GB ; second column, the edges E LB ; third column, the edges E LG ; fourth column, the combined edge maps.

the Canny algorithm can generate good edge maps, depending on the bands used. The second and third rows show the results for the MCG and HySPADE algorithms, respectively, at different threshold values in order to unveil the isoluminant edge between the granite and limestone rocks. The MCG algorithm second row) picks up the weak edge only after its tolerance is increased to a degree that results in the detection of a significant number of false edges (second row, fourth column). On the other hand, HySPADE offers a less-noisy edge map compared to the MCG algorithm; nonetheless, the background-granite and granite-limestone edges are not well defined, as shown in the third row, fourth column. Moreover, the high computational cost of the HySPADE algorithm makes it hard for the user to fine tune its tolerances, which is a clear disadvantage of the HySPADE algorithm. (More details

D ATASET C ONTAINING B, G AND L C LASSES (R AW D ATA )

regarding computational costs are given in Section IV). We also observe that at the cost of a slight increase in the number of false edges, the SRC algorithm can clearly define the background-granite edge with respect to the granite-limestone edge (fourth row, fourth column). Finally, the results of the ASRC algorithm (fifth row, fourth column) are better than all the previous algorithms in terms of clearly defining both the strong and weak edges. The ASRC algorithm also discards all of the false edges in the the background. By utilizing the available ground-truth information for the DWELL datasets, we derived reference edge maps for the scenes under study. These edge maps are utilized to compute the empirical detection and false-alarm probabilities, PD and PF , respectively, for the five algorithms (Canny applied on different bands, MCG, HySPADE, SRC and ASRC). The detection probability (also known as the sensitivity of the algorithm) corresponds to the probability that an actual edge (provided by the ground truth) is detected by the algorithm under evaluation. The false-alarm probability (also known as the complement of the specificity of the algorithm) is the probability that the algorithm detects a non-existing edge. For each algorithm, we have tuned the respective parameters in order to unveil the isoluminant edges (the assessment was made by visual inspection). We have conditioned the algorithms’ parameters to detect isoluminant edges because they present one of the most challenging problems in multicolor edge detection. The metrics PD and PF were computed by comparing the ground-truth edge-map with the algorithm outcome on a pixel-by-pixel basis. From the results presented in Table III we see that the best performance achieved by the Canny algorithm is when it is applied to band 9 (PD = 0.4533 and PF = 0.0082). It is important to note that the Canny algorithm, applied to this band, is capable to partially detect the isoluminant edge (see Fig. 5, top-right). However, without previous knowledge of the scene and the results of the application of the Canny algorithm to every band, it would be difficult to guess which band gives the best results. The MCG algorithm, on the other hand, cannot detect the isoluminant edges without producing a high number of false edges. Indeed, when the isoluminant edge is detected (second row, fourth column) the MCG performance is given by a high detection (PD = 0.9600) but also with a high falsealarm probability (PF = 0.6112). At the cost of a tremendous increase of computation complexity (see Section IV), the HySPADE algorithm outperforms the Canny algorithm in terms of sensitivity (PD = 0.7867) and the MCG algorithm

2324

IEEE TRANSACTIONS ON IMAGE PROCESSING, VOL. 23, NO. 5, MAY 2014

Fig. 6. Comparison among the Canny algorithm applied to individual bands 1, 6, 8 and 9 (first row), MCG algorithm (second row), HySPADE algorithm (third row), SRC (fourth row) and ASRC (fifth row) for the dataset containing Phyllite, Granite and Limestone rocks (second dataset). The MCG and HySPADE results are presented for a sequence of increasingly permissive tolerances in order to unveil the isoluminant edge. Last two rows show the SRC and ASRC edge maps: first column, the edges E PG ; second column, the edges E PL ; third column, the edges E LG ; fourth column, the combined edge maps. TABLE IV C OMPARISON TABLE FOR THE R ESULTS OF F IVE A LGORITHMS (C ANNY, MCG, H Y SPADE, SRC AND ASRC) FOR THE D ATASET C ONTAINING B, P, G AND L C LASSES (N ORMALIZED D ATA )

in terms of low false alarm probability (PF = 0.0565). In contrast, the SRC algorithm outperforms the previous algorithms in terms of both simplicity and sensitivity with PD = 0.9467, at the cost of a slight increase in the false-alarm probability (PF = 0.0862) in comparison to HySPADE. The ASRC algorithm outperforms all the other four algorithm in terms of highest detection and lowest false-alarm probability, PD = 0.9733 and PF = 0.0244, respectively. The edge-detection results for the second DWELL dataset are presented in Fig. 6 for intensity-normalized data. Table IV summarizes the detection and false-alarm probabilities achieved by each one of the five algorithms for this dataset. The second dataset is more challenging than the first dataset because the two classes with the isoluminant edge (i.e., granite and limestone rocks) are now positioned against a phyllite backdrop that exhibits less contrast than the blackbody. Moreover, the data is intensity normalized.

As before, the Canny edge detector achieves good performance when applied to band 9 (PD = 0.7854 and PF = 0.0301). It is very interesting to note that (for this band) the Canny algorithm is capable detecting the isoluminant edge between the granite and limestone rocks almost fully. This is because the normalization process smooths some intensity peaks and improves the contrast between granite and limestone (for this particular band) as a secondary effect. This result proves that the first category of algorithms (those that do not use spectral information) can achieve good detection as long as the best band is identified through pre-processing of the data, which can be a very difficult requirement. As for the MCG-generated edge maps, Fig. 6 (second row), the weak edge is detected only when the false-alarm probability reaches unacceptable levels. The HySPADE algorithm performs worst than the MCG algorithm (PD = 0.7445 and PF = 0.0833) and it is not capable of detecting the isoluminant edge. In contrast, the SRC algorithm recovers the strong edges as well as the weak edge between the granite and limestone rocks. Indeed, Fig. 6 (fourth row) shows a high-resolution weak edge captured by the SRC algorithm. The achieved detection and false-alarm probabilities of the SRC algorithm (PD = 0.8593 and PF = 0.0873) corroborate this observation. It is important to note that event though the SRC algorithm is able to detect isoluminant edges for challenging scenarios, it still suffers from detecting false edges for each pair of materials, as observed in both examples. However, the ASRC algorithm reduces the detection of false edges substantially (PF = 0.0652 for ASRC compared to PF = 0.0873 for SRC), owing to the fusion of material classification in the edge-detection process. The ASRC is also able to improve the detection of edges, as noted by the improved detection probability (PD = 0.8919 for ASRC compared to PD = 0.8593 for SRC). From these results, we can conclude that the SRC algorithm outperforms the MCG and HySPADE algorithms for the task of detecting edges using spectral data with minimal intensity contrast. Moreover, it performs as good as the Canny edge detector without the difficult requirement for pre-selecting the optimal band. Moreover, at the cost of a slight increase in computational cost, the ASRC algorithm outperforms all other four algorithms presented in this paper. Next, we compare the multicolor algorithms (SRC, ASRC, MCG and HySPADE) in terms of their computational costs. IV. C OMPLEXITY A NALYSIS In this section, we estimate the complexity of the feature extraction stage in the SRC and ASRC algorithms and compare it to those for the MCG and the HySPADE algorithms. Since the edge signature identification is made offline and before the edge-identification stage, we do not include its computational cost. For simplicity, in the cost estimates we will regard the cost of all operations (e.g., multiplication, addition, etc.) as equal. The SRC operations per pixel include the 2M R ratios required to form the matrix KAB (u)(i,

j ) in(4), plus the 4M R computations required to form  KAB (u) (i, j ) in (6), plus

PASKALEVA et al.: MODEL-BASED EDGE DETECTOR FOR SPECTRAL IMAGERY

TABLE V C OMPARISON TABLE FOR THE T OTAL N UMBER OF O PERATIONS R EQUIRED FOR THE SRC, ASRC, MCG AND H Y SPADE A LGORITHMS

the M R operations required to define the edges in (8). The total number of operations for the SRC algorithm is 7M R operations per pixel. The ASRC computations include those from the SRC algorithm (7M R operations per pixel) plus i, j those required to compute and utilize the parameters γm,r . To calculate these parameters, we first require the classification and label comparison of the pixels within the mask, a task that will cost 2M +6R operations per pixel. Next, the computation i, j of the γm,r parameter requires 11M R operations per pixel (two XOR operations and one OR operation for each entry in (6)). The total number of operations for the ASRC algorithm is therefore 2M + 6R + 18M R operations per pixel. Meanwhile, the MCG algorithm requires 10K −3 operations to compute the first fundamental form for each hyper-pixel, nine operations to compute the corresponding eigenvalues, and three operations to compute the monitor function and apply the threshold. The total number of operations for the MCG algorithm is therefore 10K + 9 operations per pixel. Next, for each hyper-pixel, the HySPADE algorithm requires the computation of I J spectral angles (each spectral angle costs 6K + 1 operations), plus the 2I J + 1 operations per pixel of the SA-cube to compute the one-dimensional derivative approximation, plus the I J + 1 operations required to account for the statistical accumulation of each pixel within the SA-cube. The total number of operations for the HySPADE algorithm is therefore I J (3I J + 6K + 3) operations per pixel. In Table V we present a summary of the estimated values for the four algorithms considering the same AHI and DWELL experiments we previously discussed in Section III. The proposed algorithms do not change their respective computational costs for the two examples because the edge signature identification removes the dependency of the algorithms on the actual number of bands of the data. From the presented table, we can observe that the SRC algorithm gives a 71 fold gain in computational efficiency over the MCG algorithm for the two class edge detection problem over the AHI data, whereas the ASRC algorithm gives a 23 fold gain. V. C ONCLUSIONS We have introduced two model-based, spatio-spectral edgedetection algorithms, termed the SRC and ASRC algorithms. The SRC algorithm enables the detection of edges that are due to either material change or intensity variation in scenes containing a prescribed set of materials. The ASRC algorithm

2325

is a specialized version of the SRC algorithm, aimed at detecting edges that are due to a change in the material only. The ASRC aims to reduce the detection of false edges due to unwanted changes in the intensity. Both algorithms utilize spectral library information to construct a sparse, non-separable and 3D edge operator while exploiting the concept of spectral ratio contrast. The reported SRC edge detector performs as well as the MCG edge detector for moderately challenging edges, with the advantage of requiring less operations than that required by the MCG algorithm (a reduction by a factor of 71 in our examples). However, for challenging imagery containing isoluminant edges, the SRC and ASRC edge detectors outperform the MCG and HySPADE edge detectors by a wide margin, as quantified by the detection and false-alarm probabilities. This provides a strong validation of the efficacy of the spectral ratio contrast concept by showing that the use of select band ratios can lead to reliable identification of weak edges in the presence of noise. Moreover, with a slight increase in the complexity (3 folds in our examples), the ASRC algorithm, which also involves classification-based step, is capable of minimizing the false-alarm edges, outperforming the SRC, MCG and HySPADE algorithms. The dramatic reduction in the number of operations with respect to other algorithms such as the MCG and the HySPADE algorithms is a key advantage of the proposed algorithms. The reduced number of operations is mainly due to the property that only a few bands are required to perform edge detection. In principle, this property can be exploited to reduce the spectral-image acquisition time substantially by requiring the sensing of only those bands that are most relevant to the set of materials within the scene. In particular, this band-reduction feature is particularly relevant to emerging spectral imaging sensors that are bias tunable, such as the DWELL sensor, where one can perform intelligent acquisition by programming the sensor electronically to sense only at the few prescribed bands. R EFERENCES [1] R. Gonzalez and R. Woods, Digital Image Processing. Upper Saddle River, NJ, USA: Prentice-Hall, 2001. [2] J. Canny, “A computational approach to edge detection,” IEEE Trans. Pattern Anal. Mach. Intell., vol. 8, no. 6, pp. 679–698, Nov. 1986. [3] I. Sobel and G. Feldman, “A 3×3 isotropic gradient operator for image processing,” Ph.D. dissertation, Stanford Artif. Project, Stanford Univ., Stanford, CA, USA, 1968. [4] P. Blomgren and T. F. Chan, “Color TV: Total variation methods for restoration of vector valued images,” IEEE Trans. Image Process., vol. 7, no. 3, pp. 304–309, Mar. 1998. [5] L. Lucchese and S. K. Mitra, “Color image segmentation: A state-ofthe-art survey,” Proc. Ind. Nat. Sci. Acad., vol. 67, no. 2, pp. 207–221, 2001. [6] X. Chen and H. Chen, “A novel color edge detection algorithm in RGB color space,” in Proc. IEEE 10th ICSP, Oct. 2010, pp. 793–796. [7] G. Sapiro, “Color snakes,” Comput. Vis. Image Understand., vol. 68, no. 2, pp. 247–253, 1997. [8] B. Sandberg and T. F. Chan, “Logic operators for active contours on multi-channel images,” Dept. Math., UCLA, Los Angeles, CA, USA, CAM Rep. (02-12), 2002, pp. 2–12. [9] P. Bao, D. Zhang, and X. Wu, “Canny edge detection enhancement by scale multiplication,” IEEE Trans. Pattern Anal. Mach. Intell., vol. 27, no. 9, pp. 1485–1490, Sep. 2005.

2326

[10] W. Jiang, K.-M. Lam, and T.-Z. Shen, “Efficient edge detection using simplified Gabor wavelets,” IEEE Trans. Syst., Man, Cybern. B, Cybern., vol. 39, no. 4, pp. 1036–1047, Aug. 2009. [11] Y. Jing, J. An, and Z. Liu, “A novel edge detection algorithm based on global minimization active contour model for oil slick infrared aerial image,” IEEE Trans. Geosci. Remote Sens., vol. 49, no. 6, pp. 2005–2013, Jun. 2011. [12] P. Xu, Q. Miao, C. Shi, J. Zhang, and M. Yang, “General method for edge detection based on the shear transform,” IET Image Process., vol. 6, no. 7, pp. 839–853, Oct. 2012. [13] S. Di Zenzo, “A note on the gradient of a multi-image,” Comput. Vis. Graph. Image Process., vol. 33, no. 1, pp. 116–125, 1986. [14] A. N. Evans and X.-U. Liu, “A morphological gradient approach to color edge detection,” IEEE Trans. Image Process., vol. 15, no. 6, pp. 1454–1463, Jun. 2006. [15] J.-F. Rivest, P. Soille, and S. Beucher, “Morphological gradients,” J. Electron. Imag., vol. 2, no. 4, pp. 326–336, 1993. [16] P. Soille, Morphological Image Analysis: Principles and Applications. Berlin, Germany: Springer-Verlag, 2003. [17] A. Cumani, “Edge detection in multispectral images,” CVGIP: Graph. Models Image Process., vol. 53, no. 1, pp. 40–51, 1991. [18] G. Sapiro and D. Ringach, “Anisotropic diffusion of multivalued images with applications to color filtering,” IEEE Trans. Image Process., vol. 5, no. 11, pp. 1582–1586, Nov. 1996. [19] P. E. Trahanias and A. N. Venetsanopoulos, “Color edge detection using vector order statistics,” IEEE Trans. Image Process., vol. 2, no. 2, pp. 259–264, Apr. 1993. [20] P. J. Toivanen, J. Ansamäki, J. P. S. Parkkinen, and J. Mielikäinen, “Edge detection in multispectral images using the self-organizing map,” Pattern Recognit. Lett., vol. 24, no. 16, pp. 2987–2994, 2003. [21] J. Jordan and E. Angelopoulou, “Edge detection in multispectral images using the n-dimensional self-organizing map,” in Proc. 18th IEEE ICIP, Sep. 2011, pp. 3181–3184. [22] M. Yu, G. Xiaodong, and W. Yuanyuan, “Feature fusion method for edge detection of color images,” J. Syst. Eng. Electron., vol. 20, no. 2, pp. 394–399, 2009. [23] A. Mittal, S. Sofat, and E. Hancock, “Detection of edges in color images: A review and evaluative comparison of state-of-the-art techniques,” in Autonomous and Intelligent Systems (Lecture Notes in Computer Science), M. Kamel, F. Karray, and H. Hagras, Eds. Berlin, Germany: Springer-Verlag, 2012, pp. 250–259. [24] A. Koschan and M. Abidi, “Detection and classification of edges in color images,” IEEE Signal Process. Mag., vol. 22, no. 1, pp. 64–73, Jan. 2005. [25] E. Nezhadarya and R. K. Ward, “A new scheme for robust gradient vector estimation in color images,” IEEE Trans. Image Process., vol. 20, no. 8, pp. 2211–2220, Aug. 2011. [26] R. G. Resmini, “Hyperspectral/spatial detection of edges (hyspade): An algorithm for spatial and spectral analysis of hyperspectral information,” Proc. SPIE, vol. 5425, pp. 433–442, Jan. 2004. [27] J. Weier and D. Herring. (2010). Measuring Vegetation (NDVI and EVI) [Online]. Available: http://earthobservatory.nasa.gov/Features/Measuring Vegetation/ [28] A. Huete, “A soil-adjusted vegetation index (SAVI),” Remote Sens. Environ., vol. 25, no. 3, pp. 53–70, 1988. [29] Y. Kaufman and D. Tanre, “Atmospherically resistant vegetation index (ARVI) for EOS-MODIS,” in Proc. IEEE Int. Geosci. Remote Sens. Symp., Mar. 1992, pp. 261–270. [30] K. Mayeda, A. Hofstetter, J. O’Boyle, and W. Walter, “Stable and transportable regional magnitudes based on coda-derived moment-rate spectra,” Bull. Seismol. Soc. Amer., vol. 93, pp. 224–239, Feb. 2003. [31] K. Mayeda, L. Malagnini, and W. Walter, “A new spectral ratio method using narrow band coda envelopes: Evidence for non-self-similarity in the Hector mine sequence,” Geophys. Res. Lett., vol. 34, no. 11, p. L11303, 2007. [32] H. Hartse, S. Taylor, W. Phillips, and G. Randall, “Regional seismic discrimination in central Asia with emphasis on Western China,” Los Alamos Nat. Lab., Los Alamos, NM, USA, Tech. Rep. LA-UR 96-2390, 1996. [33] H. J. Patton and B. J. Mitchell, Monitoring the Comprehensive NuclearTest-Ban Treaty: Regional Wave Propagation and Crustal Structure. Cambridge, MA, USA: Birkhäuser-Verlag, 2001. [34] U. Al-Suwailem and J. Keller, “Multichannel image identification and restoration using continuous spatial domain modeling,” in Proc. Int. Conf. Image Process., vol. 2. 1997, p. 466.

IEEE TRANSACTIONS ON IMAGE PROCESSING, VOL. 23, NO. 5, MAY 2014

[35] P. G. Lucey, T. J. Williams, J. L. Hinrichs, M. E. Winter, D. Steutel, and E. M. Winter, “Three years of operation of AHI: The University of Hawaii’s Airborne Hyperspectral Imager,” Proc. SPIE, vol. 4369, pp. 112–120, Oct. 2001. [36] J. Andrews et al., “Demonstration of a bias tunable quantum dots-in-awell focal plane array,” Infr. Phys. Technol., vol. 1, pp. 1–5, Jun. 2009. [37] M. Loog and B. Ginneken, “Segmentation of the posterior ribs in chest radiographs using iterated contextual pixel classification,” IEEE Trans. Med. Imag., vol. 25, no. 5, pp. 602–611, May 2006. [38] T. Rohlfing and C. R. Maurer, “Multi-classifier framework for atlasbased image segmentation,” Pattern Recognit. Lett., vol. 26, no. 13, pp. 2070–2079, 2005. [39] B. Pasakaleva, W.-Y. Jang, M. Hayat, Y. Sharma, S. Bender, and S. Krishna, “Multispectral classification with bias-tunable quantum dots-in-a-well focal plane arrays,” IEEE Sensors J., vol. 11, no. 6, pp. 1342–1351, Jun. 2011. [40] D. A. B. Miller et al., “Band-edge electroabsorption in quantum well structures: The quantum-confined Stark effect,” Phys. Rev. Lett., vol. 53, pp. 2173–2176, Nov. 1984.

Biliana S. Paskaleva received the B.S. degree in electrical engineering from the Technical University, Varna, Bulgaria, the M.S. degree in electrical engineering from the Department of Electrical and Computer Engineering, University of New Mexico, Albuquerque, and the Ph.D. degree in electrical engineering from the University of New Mexico, under the supervision of Prof. M. M. Hayat, in 1992, 2004, and 2010, respectively. She is a Senior Member of Technical Staff with the Sandia National Laboratories, Albuquerque. Her research interests are in the area of remote sensing and image processing, spectral and spatial feature extraction, multimodal information retrieval, Bayesian statistics and neural networks, validation and verification, and uncertainty quantification.

Sebastián E. Godoy was born in Concepcin, Chile. He received the B.Sc. and M.Sc. degrees in electrical engineering from the University of Concepcin, Concepcion, in 2006 and 2009, respectively. He is currently pursuing the Ph.D. degree in electrical engineering with the Center for High Technology Materials, University of New Mexico, Albuquerque. His research interests are in the area of algorithms for recognition and classification, image processing, and feature extraction.

Woo-Yong Jang received the M.S. degree in electrical engineering from the University of Southern California, Los Angeles, and the Ph.D. degree in electrical and computer engineering from the University of New Mexico, Albuquerque, in 2004 and 2012, respectively. He was a Research Support Engineer at Solid State Scientific Corporation from 2012 to 2013. He currently holds the National Research Council Fellowship as a Research Associate with the U.S. Air Force Research Laboratory. He has over 40 peer-reviewed journals and conferences and holds one U.S. patent. His research areas of interest include compressive spectral sensing, plasmonic focal-plane array imaging, multimodal feature extraction and integration for target identification and tracking, and system engineering for sensing algorithms.

Steven C. Bender was with the Space and Remote Sensing Sciences Group, Los Alamos National Laboratory, Los Alamos, NM, USA. His research areas include remote sensing, sensor modeling, and system integration and performance.

PASKALEVA et al.: MODEL-BASED EDGE DETECTOR FOR SPECTRAL IMAGERY

Sanjay Krishna (F’14) received the master’s degree in physics from IIT Madras, Chennai, India, and the M.S. degree in electrical engineering and the Ph.D. degree in applied physics from the University of Michigan, Ann Arbor, in 1996, 1999, and 2001, respectively. He joined the University of New Mexico, Albuquerque, as a Tenure Track Faculty Member, in 2001, where he is currently a Professor with the Electrical and Computer Engineering Department, Center for High Technology Materials. He has authored and co-authored more than 200 peer-reviewed papers in journals (h-index = 29) and conferences, and holds six patents. His current research interests include growth, fabrication, and characterization of self-assembled quantum dots and type-II InAs- and InGaSb-based strained layer superlattices for midinfrared detectors. He was a recipient of various awards, including the Gold Medal from IIT Madras, the Early Career Achievement Award from SPIE, NAMBE, and the IEEE Nanotechnology Council, and the UNM Teacher of the Year Award. He is a fellow of IEEE and SPIE.

2327

Majeed M. Hayat was born in Kuwait in 1963. He received the B.Sc. (summa cum laude) degree in electrical engineering from the University of the Pacific, Stockton, CA, USA, and the M.S. and Ph.D. degrees in electrical and computer engineering from the University of Wisconsin-Madison, in 1985, 1988, and 1992, respectively. He is currently a Professor of Electrical and Computer Engineering, an Associate Director of the Center for High Technology Materials, and a General Chair of the Optical Science and Engineering Program with the University of New Mexico. His research activities cover a broad range of topics, including avalanche photodiodes, signal and image processing, algorithms for spectral and radar sensing and imaging, optical communication, networked computing, and modeling interdependent networks with applications to smart grids. He was an Associate Editor of Optics Express from 2004 to 2009 and the Chair of the topical committee of Photodetctors, Sensors, Systems and Imaging with the IEEE Photonics Society from 2010 to 2013. He was a recipient of the National Science Foundation Early Faculty Career Award in 1998 and the Chief Scientist Award for Excellence in 2006 from the National Consortium for MASINT Research of the Defense Intelligence Agency. He has authored and co-authored over 85 peer-reviewed journal articles (over 2800 citations, h-index = 28), and has eight issued patents, three of which have been licensed. He is a fellow of IEEE, SPIE, OSA.