In this article we propose an algorithm for multichannel image segmentation within the prede ned region of interest. The algorithm allows to determine whether the region of interest is homogeneous or not and it also produces a partition of the region of interest into homogeneous subsets and re nes objects borders. The proposed method has been applied to remote sensing images of agricultural elds as a preliminary step of crop detection mechanism. Segmentation quality has been estimated by means of a speci c measure proposed in this article.
There are a lot of applications that require automatic segmentation to be
performed as initial step [
Existing image segmentation methods can be regarded in the scope of
further classi cation of popular segmentation techniques given in [
In our article we assume that the database of objects' contours is available. This database may be presented in the form of vector map from geoinformation system (GIS). We explore the case of automatic analysis of internal map object structure and its border re nement, thus we consider only the explored object and its surroundings as input data. The aim of the segmentation in our case is to de ne signi cant object's parts that di er from each other by their spectral-spatial characteristics, i.e. to clarify an internal object structure. Resulting contours should separate an initial object into meaningful parts. The outside object border is supposed to be known from the vector map and is referred as the region of interest (ROI). It is assumed that multiple features are available.
As far as we suppose that it is more important to identify borders inside the
ROI rather than to re ne its contour, we have based our algorithm on a merge
approach [
Managing the splitting and merging process allows our algorithm to work adaptively with each particular contrast and texture, which is highly determined by image resolution. Such opportunity is very important for the agricultural elds, because usually images are very di erent not only in sense of planted crops but also as a result of relief and image acquisition properties. Relief changes lead to texture variations within the elds with the same crop class. And acquisition parameters signi cantly determine particular image contrast properties.
The proposed algorithm has been revised by means of segmentation quality measure on the set of remote sensing images. Experiments have shown successful results in most of the cases and the best pair of the algorithm parameters has been recommended.
The article includes three main parts: problem statement, algorithm description and experimental evaluation of the method on remote sensing images. The article ends with conclusions and acknowledgements sections. 2
In our notations initial image of L features with size N1 symbol Xl (n1; n2) ; 0 n1 N1 1; 0 n2 N2 1; 0 and n2 are the coordinates of image pixel. The general segmentation problem is formulated as optimization problem:
S = fSkg ; 0 k
KS 1; = f ig ; 0 i
K 1; [ Sk ;
k Q (S; ) 7! max: (1) (2) where Q (S) is the segmentation quality measure, S = fSkg ; 0 k KS 1 is a set of pixel subsets de ning the results of segmentation, KS is the number of subsets, = f ig ; 0 i K 1 is an ideal segmentation and K is ideal number of subsets.
We o er the following measure segmentation quality evaluation:
Q (S; ) =
where j:::j is the number of pixels per subset. The segmentation quality measure
is based on the quantities introduced in [
Let us denote the source remote sensing image as Fm (n1; n2) ; 0 m M 1; 0 n1 N1 1; 0 n2 N2 1. Here by we describe all of the processing steps of our algorithm including preprocessing and feature calculation.
Step 1. Input data preparation.
We clip initial remote sensing image within the bounding box of vector map object to be analyzed using standard GIS instruments. The resulting piece of the image represents the region of interest and contains all pixels corresponding to analyzed object and its surroundings.
Step 2. Feature extraction.
This stage includes calculation of initial feature image Xl (n1; n2) ; 0 n1 N1 1; 0 n2 N2 1; 0 l L 1 using local averages and variations for each pixel block with W W size in the source image. The number of features L = 2M is two times higher than the original image channel count, because these features are computed separately per each image component.
Step 3. Feature variance normalization.
To exclude the impact of di erent dynamic ranges for the computed features we normalize each component using its global standard deviation l: xl (n1; n2) =
Xl (n1; n2) : l (4) Step 4. Splitting image into tiny regions.
4.1. Before splitting initial image with segment indexes I (n1; n2) ; 0 n1 N1; 0 n2 N2 is initialized by zero: I (n1; n2). Current number of segments is set to zero S = 0 too. Further analysis is produced for the sequence of grids de ned by shifting parameter a = 1; 4.2. For each feature vector in the grid x (n1; n2) ; (n1; n2) 2 T with step aW : T = f(n1; n2) : n1 mod aW = 0; n2 mod aW = 0g ; (5) where mod is the remainder after division, the neighborhood is determined. This neighborhood contains pixels located at distance aW per each spatial coordinate from the current pixel. There are three possible cases. First case is when the neighborhood pixels have di erent indexes of segments. Second case is when only one pixel has non zero segment index s for the current neighborhood. And the last case is when aW = 1 and all neighborhood pixels have the same index s. If the second and third situations take place, then the condition of the similarity between feature vectors is checked: (x (n1; n2) ; x (n1 aW; n2 aW )) < "; (n1; n2) 2 T (6) where is some distance metric i.e. Euclidean distance. If the condition (6) is ful lled, the segment index I (n1; n2) is set to s. Otherwise, it is assigned to smax + 1, where smax is the biggest segment index at the current time.
In the rst case, when there are several di erent values of segment index in the neighbourhood, current feature vector have to be compared with an average feature vectors for each of the presented segments: (7) (8) (9) (x (n1; n2) ; E (s)) < "; (n1; n2) 2 T; s = argmins2fI(n1 aW;n2 aW )g (x (n1; n2) ; E (s)) where E (s) is an average feature vector for the segment with index s. If this condition is ful lled for some index s, it is stored in I (n1; n2) as a current estimate of the segment number.
4.3 The grid step is decreased in twice a = a2 and the steps from 4.1 to 4.3 are repeated while the condition aW 1 is valid.
Step 5. Merging.
For each pair of the adjacent segments two segments are merged if their average feature vectors are closer than the threshold value:
(E (s1) ; E (s2)) < " Step 6. ROI border re nement.
This step is intended to the analysis of non zero points on the segment border. This analysis includes three steps listed below:
6.1 Distance map calculation for each feature vector. This map contains the distances between current feature vector and average feature vectors of the segments corresponding to pixels within the local window W W .
6.2 Parabolic ltering of the distance map. Each value of the distance within the window is multiplied by the coe cients of parabolic lter with same size W W and then the minimum is used as a result of ltration for current window position. The coe cients of parabolic lter are de ned by following equation:
P (w1; w2) = w1 w1 + w2 w2; (10) where w1 and w2 take values from the range W2 ; W2 .
6.3 Final index of the segment is de ned as the index of segment with the minimum ltered distance which is also less than ".
The algorithm produces an image with the segments' indexes. Step 2 can be modi ed if it is necessary to apply any other features. 4
The proposed algorithm has been used for agricultural elds analysis. Input images were acquired by sensor UK DMC for the Samara region. Image resolution was 22 meters per pixel. Test images were obtained during the period from the 1st of April 2012 till the 15th of June.
The aim of our experimental research was to de ne which parameter values deliver the best method performance. The following groups of factors may signi cantly a ect segmentation quality. The rst group accumulates images which are partially clouded and with shadows from the clouds. For these images clouds and shadows always are distinguished as di erent segments and leads to over segmentation. This is irremovable error, therefore only cloudless images should be used. The second group accumulated images of elds with the rough relief structure. If there are some ditches, channels, ravines or the slope and curvature is high, it would be di cult to make segmentation properly. The last group includes rough-textured elds. Usually, non hybrid eld has the rough texture of its surface because of the agricultural management events. In this case there can be regular lines or spots on the images inside the borders of the eld. If the texture elements are enough big, they may lead to over segmentation too. Thus, the results of the algorithm are highly dependent on input data and can not be applied without calibration.
For the agricultural eld border monitoring we have used infrared, red and green spectral channels. The features were local averages and variances for each spectral channel.
We have used a sample set of 120 images for the typical problem cases listed above (30 images per category) and 30 images with simple non hybrid structure (as an example of the case of good conditions for the algorithm). The ideal segmentation has been done manually. The examples of test images are presented in g. 1.
The quality of segmentation for each image has been estimated by means of Q criterion described above. We have tested algorithm parameters for the a) b) c) d)
To inspect the relationship between each category and quality of segmentation we have assigned the index g to each of the four image categories: g = 1 for hybrid eld, g = 2 for non hybrid eld with complex relief, g = 3 for simple non hybrid eld, g = 4 for textured hybrid eld. Image index in the particular group is denoted as j and there are J = 30 images per each category.
For all pairs of parameters W and " an averages and variances of Q per group and for the whole image set have been calculated using formulas: mQg (W; ") =
Qgj (W; ") ;
J 1 X J
j=1 v
J Qg (W; ") = tuu J1 X (Qgj (W; ") j=1
mQg (W; "))2; mQ (W; ") =
For each category we have selected an optimal pair of parameters Poptg as a pair that ful ll the following condition:
Poptg = f(W; ") : mQg (W; ") mg gg ; where mg = max(W;") [mQg (W; ")], g = max(W;")2Pg [ Qg (W; ")], (11) (12) (13) (14) (15) Pg = f(W; ") : mQg (W; ") = mgg.
The table of mg and g values for each category is shown below. Category g hybrid eld 1 non hybrid eld with complex relief 2 simple non hybrid eld 3 textured hybrid eld 4
The optimal pairs of the parameters for the whole image set have been found as intersection of optimal parameters sets for each category Poptg .
The table 2 contains optimal for the whole sample set parameters and average segmentation quality measure mQ (W; ") for them. mQ1 (W; ") is an average value of the segmentation quality measure only for the hybrid eld samples. Fig. 2. Segmentation results with the best pair of parameters W = 5 and " = 0:6 for the images from the g. 1
In this article, we have proposed an algorithm that helps to make a segmentation of particular objects. It is based on two steps: excess segmentation and merging. The quality of segmentation depends on the similarity threshold and initial degree of over segmentation that is controlled by the window size. Algorithm has been tested on agricultural eld images and the best pair of parameters have been selected. The control parameters of the algorithm allow to use it for another kind of data with di erent resolution and contrast properties, but the calibration must be carried out in advance.
Acknowledgements This study was nancially supported by RFBR projects 16-37-00043 mola, 16-29-09494 o m.