Postharvest Biology and Technology 128 (2017) 33–43
Contents lists available at ScienceDirect
Postharvest Biology and Technology journal homepage: www.elsevier.com/locate/postharvbio
Multisensor X-ray inspection of internal defects in horticultural products M. van Daela,* , P. Verbovena , J. Dhaeneb , L. Van Hoorebekeb , J. Sijbersc , B. Nicolaia a b c
Flanders Centre of Postharvest Technology/BIOSYST-MeBioS, KU Leuven—University of Leuven, Willem de Croylaan 42, 3001 Leuven, Belgium UGCT, Dept. Physics and Astronomy, Ghent University, Proeftuinstraat 86, 9000 Ghent, Belgium iMinds-Vision Lab., Department of Physics, University of Antwerp, Universiteitsplein 1, B-2610 Wilrijk, Belgium
A R T I C L E I N F O
Article history: Received 13 October 2016 Received in revised form 7 February 2017 Accepted 7 February 2017 Keywords: Non-destructive Inline Horticultural product X-ray Shape model
A B S T R A C T
A combination of 3-D vision and X-ray radiography is proposed to enable low-cost, generally applicable online inspection of internal quality of horticultural and potentially other products. The underlying concept assumes that the shape of the product is known beforehand through a deformable shape model. A 3-D vision system is used in combination with the shape model to accurately determine the complete outer surface shape of the sample. This shape is voxelized to generate a reference product from which a X-ray radiograph is simulated to be compared with a measured radiograph, hence revealing the presence of any defects or disorders. Advantages of this method are that small deviations in internal density are detected easily since the cumulative information of the bulk object shape is removed. Furthermore, no specific detection algorithms have to be developed for different types of defects, since the method will directly identify deviations from the ideal. Validation on two datasets and comparison with two reference detections methods (classical image processing and a human operator) shows that the proposed method reliably (accuracy >99%) detect defects larger than 3.5 mm radius with densities differences between sample and defects as small as 10%. Voids are reliably (accuracy >99%) detected down to a radius of 1.5 mm, corresponding to a volume of less than 0.02 cm3. © 2017 Elsevier B.V. All rights reserved.
1. Introduction Horticultural products often develop internal defects, cracks or cavities that reduce their commercial value. Currently the most widespread method to detect these is destructive evaluation where a number of samples is randomly selected from each batch and cut open for inspection (Shewfelt et al., 2012). This has some inherent disadvantages like financial losses because of destroying a number of samples and the fact that only a small subset of all products is checked. Recent research shows that a number of internal defects can be detected using non-destructive technologies such as visible/nearinfrared (NIR) spectroscopy (Bobelyn et al., 2010; Magwaza et al., 2012; Nicolaï et al., 2007), nuclear magnetic resonance (MRI, Nicolaï et al., 2014; Zhang and McCarthy, 2013), hyperspectral imaging (Haff et al., 2013; López-Maestresalas et al., 2016), and magnetic resonance imaging (MRI, Clark et al., 1999; Herremans et al., 2014; Lammertyn et al., 2003a, 2003b). X-ray based methods,
* Corresponding author. E-mail address:
[email protected] (M. van Dael). http://dx.doi.org/10.1016/j.postharvbio.2017.02.002 0925-5214/© 2017 Elsevier B.V. All rights reserved.
however, have as advantage their sensitivity to spatial density differences inside an object and the excellent penetration properties of X-rays in horticultural products (Kotwaliwale et al., 2014). X-rays populate the electromagnetic spectrum between 10 and 0.01 nm. When passing through an object, these X-rays physically interact with the object material and are partly absorbed and scattered (Barrie Smith and Webb, 2010). The remaining X-rays are recorded on a detector, resulting in an X-ray radiograph, an image containing superimposed information or a projection of the 3-D object in a 2D plane. X-ray radiography is a fast and cheap 2D imaging technique and is already used in industry for detection of foreign bodies with distinct density (e.g. metals, stones) such that they appear with high contrast in the resulting image for easy segmentation and detection. Internal defects usually cause smaller attenuation differences in the object. Furthermore, radiographs render the cumulative attenuation along the X-ray path from source to detector and therefore are affected by the size and shape of the scanned object. Contrast gradients in the image are thus not only due to density variations due to defects and other local features, but are also caused by object shape and size. Although radiographs have proven to be useful in detecting certain types of
34
M. van Dael et al. / Postharvest Biology and Technology 128 (2017) 33–43
internal defects in horticultural products (Haff et al., 2006; Hansen et al., 2005; van Dael et al., 2016a), variations in their shape, size, position and density reduce the classification accuracy of the method, especially when defects are small. Shape variation can be considerable for agricultural and horticultural products (Rogge et al., 2015). In X-ray computed tomography (CT) multiple radiographs of the same object, taken from different angles, are combined to produce a 3-dimensional image using a mathematical algorithm (Magwaza and Opara, 2014; Nicolaï et al., 2014). The resolution varies from 1 mm in medical scanners down to the nanometer level on dedicated micro-CT scanners (Verboven et al., 2008). X-ray CT has been used successfully to detect internal disorders in pear (Lammertyn et al., 2003a), apple (Herremans et al., 2014; Nicolaï et al., 2014; Schatzki et al., 1996) and pineapple (Haff et al., 2006). However, the equipment is expensive and the 3-dimensional reconstruction comes at a high computational cost making it difficult to apply into existing sorting lines. The aim of this paper is to improve the detection accuracy of internal defects of horticultural products by means of inline X-ray radiography. To achieve this objective, the method patented in (van Dael et al., 2016b) is applied to internal quality inspection of horticultural products, explicitly taking into account the shape, size and orientation of the object by means of a 3-D vision system and a deformable shape model. This normalizes the measured radiograph for the bulk shape and size of the inspected sample by simulating a reference radiograph from the shape model that is fitted to the point cloud produced by the 3D-vision system. A range of product shapes and defects are considered. The principle of the method is demonstrated by comparing it to radiography detection algorithms which use a conventional method of imaging and detection. 2. Materials and methods 2.1. Products The shape of the majority of horticultural products can be approximated with curvilinear shapes. In this work, we have first considered a family of curvilinear shaped model objects (tori) with artificially induced defects of different sizes and contrast and subsequently a real horticultural product, namely ‘Conference’
pears with a defect known as core breakdown (Franck et al., 2007; Lammertyn et al., 2003a,b). All simulations and calculations were performed on a on a desktop pc with a quad-core 3.4 GHz processor, 16 GB of RAM memory and an nVidia Quadro K600 1GB graphics card. 2.1.1. Torus shapes A set of random tori (1 mm < r < 25 mm, 1 mm < R < 25 mm, N = 2240, with r and R the minor and major radius respectively and N the number of samples) was voxelized according to Patil and Ravi (2005) to a binary 3-dimensional matrix. These dimensions were chosen arbitrarily to not exceed a maximum diameter of 100 mm, approximating the size range of many horticultural products. Spherical defects (0.5 mm < r < 7 mm) were added at random locations with voxel values (representing fractional density of the defect to the object material) ranging from 0 (voids) to 1 (no defect present) in 0.1 increments. Fig. 1 shows 3 tori with added defect shapes in the top row, together with a radiograph simulated from a random orientation in the bottom row. 2.1.2. X-ray CT imaging of ‘Conference’ pears The realistic dataset consisted of a set of X-ray CT scans of Pyrus communis cv. ‘Conference’ pears harvested from orchards in the province of Limburg, Belgium between September 17–23, 2013. The pears were stored at 1 C and environmental gas conditions until January 7, 2014 when the X-ray CT scans were performed. Samples were randomly selected (N = 30) and scanned in a Philips AEA Tomohawk X-ray CT system equipped with a Nikon metrology 160 Xi Gun set and a Thomson TH 9428HX – Adimec MX12P image intensifier – CCD camera pair. Source parameters were 65 kV at 500 mA, 380 projections were captured with a rotation step of 0.5 and reconstructed resulting in an effective voxel size of 131 mm (Lammertyn et al., 2003a,b). A subset of the dataset is shown in Fig. 2. Radiographs from 15 random orientations were simulated from each pear (N = 450) as described in Step 2. X-ray simulation. Defects were manually segmented from CT scans and photographs were taken to serve as a ground truth. Since the segmented defects are elongated in shape along the core, they were approximated with ellipsoids using a normalized second central approach (Mukundan and Ramakrishnan, 1998). A histogram of the ellipsoid axis lengths is shown in Fig. 3 to quantify the size of the defects.
Fig. 1. Top row: three tori with added defects rendered in red. Bottom row: radiographs simulated from the randomly oriented tori with a detector pixel size of 0.5 mm using the ASTRA toolbox. Defect densities in these examples are 0% of sample densities i.e. voids.
M. van Dael et al. / Postharvest Biology and Technology 128 (2017) 33–43
35
Fig. 2. The pear dataset used for experimental validation. The first column shows a radiograph as captured by the CT system. The second and third column show a reconstructed CT slice and visual image respectively through approximately identical planes. The fourth and fifth column show simulated radiographs from random orientations. Samples in the first and second row contain no defects while those in the third and fourth row do.
2.2. Multisensor detection system A combination of 3-D vision and X-ray radiography is proposed for improved detection of internal disorders. In what follows the ‘sample’ is defined as a voxelized torus with added defect or a pear CT scan, respectively, as described above. In principle the samples can have a random orientation with respect to the scan setup. As this research reported here is a proof of concept, a virtual multisensor detection system will be considered, and the outputs of the 3-D vision system and X-ray system are obtained through simulation using the datasets described in Section 3.1:
2. The 3-D vision system generates a subset of sample surface coordinates. This is done by random selection of surface points and by adding zero-centered Gaussian noise (m = 0, s = 0.1 mm). 3. The X-ray system generates a radiograph of the sample. This measured radiograph is obtained in this study by simulation using the dedicated ASTRA toolbox (Palenstijn et al., 2013, 2011; van Aarle et al., 2015). Gaussian noise (m = 0, s = 2.5% of maximum grey level) was added to simulate detector noise. This corresponds to a signal-to-noise ratio of 40, which is much lower than that of any commercially available detector. The simulation is further explained in section Step 2. X-ray simulation.
1. The sample is randomly rotated with respect to the coordination system of the setup. 2.3. Defect detection To detect internal defects, the radiograph of a sample that is captured by the X-ray detector is compared with a radiograph that is simulated from a reference object free of defects. The reference object is determined by fitting the deformable shape model (see Step 1. Shape model fitting) to the surface point cloud obtained by a 3-D vision system with the same density distribution of the horticultural product that is to be inspected.
Fig. 3. Axis lengths of the ellipsoids with the same normalized second central moments as manually segmented defects present in the pear dataset.
2.3.1. Step 1. Shape model fitting The proposed method requires the surface shapes of the population of samples to be known in the form of a deformable shape model (Schneider and Eisert, 2009). During inspection this deformable shape model is fitted to the point cloud measured by a 3-D vision system. To simplify the fitting process, a principal component analysis (PCA) based shape model is used. The model is constructed a priori for the population of horticultural products (such as the pears under study here) subject to inspection by a manual alignment, determining the outer shape, finding corresponding points by interpolating to a grid in a spherical coordinate system and applying a PCA to their coordinates for a number of product items. This results in a mean position for every point,
36
M. van Dael et al. / Postharvest Biology and Technology 128 (2017) 33–43
representing the mean shape, and displacement vectors for every point represented by the principal, or shape, components (Danckaers et al., 2015). Any shape originating from the same population as from which the deformable model was built can be represented as a linear combination of the mean shape and the shape components as shown in Eq. (1). C ¼ SW þ m
ð1Þ
In Eq. (1) C is containing the 3-dimensional coordinates of the n points describing the shape. S is the n-by-s matrix containing the shape components or displacement from the mean vectors with s being the number of components used to create the shape model. The weight of each components is contained in w. The matrix m contains the 3-dimensional coordinates of the n points describing the mean model shape. The artificial dataset consisted of computer-generated torus shapes of which a set (r: m = 12.5 s = 3, R: m = 12.5 s = 3, N = 100) was used to generate a deformable shape model shown in Fig. 4. Two components describe all of the possible variation for a set of tori since they are defined by two parameters: their major and minor radii. When applied to Eq. (1), matrix S will be of size n-by-2 since only 2 shape components are used, with n the number of points describing the tori. Instances of the resulting shape model with varying weights (s, 0, s) for the first two components are shown in Fig. 4. For the pears dataset the corresponding points to build the deformable shape model were extracted from a subset of these scans by combining a manual alignment and resampling of the surface voxels coordinates to a coarse grid of spherical coordinates with origin in the pear centroid. Instances of the resulting shape model with varying weights (1 s, 0, 1 s) for the first two components are shown in Fig. 5. The shape model is fitted to the measured surface point cloud obtained by the (virtual) 3-D vision system. Various methods have been described in literature to register deformable shape models to unstructured point clouds. In this work an iterative closest point (ICP) registration (Besl and McKay, 1992) is followed by a least squares based fit of each model vertex to the closest measured point. A least squares approach is possible in the case of PCA-based models where any shape can be described as a linear combination of the mean shape and the shape components (Eq. (1)). The system is solved iteratively until convergence is reached, resulting in a transformation matrix describing the shape model’s pose and a matrix w containing the shape components weights, as described
Fig. 4. Instances of the torus shape model obtained by different linear combinations of the mean shape with the first two components. The weights of the components are 1, 0 and 1 standard deviations from left to right and top to bottom for the first and second component respectively. The first component describes the minor radius, while the second component describes the major radius. The color scale indicates the deformation from the mean shape in mm.
Fig. 5. Instances of the pear shape model obtained by different linear combinations of the mean shape with the first two components. The weights of the components are 1, 0 and 1 standard deviations from left to right and top to bottom for the first and second component respectively. The first component describes mostly sample size, while the second components describes vertical elongation. The color scale indicates the deformation from the mean shape in mm.
in Eq. (1). Fig. 6 shows the pear shape model fitted against three different point clouds and bar charts displaying the weights of the first 15 model components explaining 99% of variance. 2.3.2. Step 2. X-ray simulation The ASTRA toolbox was used to generate a simulated radiograph without defects of the fitted shape obtained in step 1. To achieve a radiograph that can be compared with the measured one, the following aspects need to be taken into account. First, the fitted shape needs to be converted from a surface mesh to a volumetric representation, a process called voxelization as shown in Fig. 7. This volumetric representation is equivalent to a 3-dimensional image, where each voxel value (a voxel is the 3-dimensional equivalent of a pixel) represents the local sample density which is in this case proportional to the local X-ray linear attenuation coefficient, The latter is the actual output of a full CTscan. Voxelization of the shape model, which essentially is a set of point coordinates and the facets connecting them, is carried out by a ray intersection algorithm as described by Patil and Ravi (2005). It passes rays along the X-axis in a preset order (incremented along the Y-axis first, then along the Z-axis). If the shape model is errorfree, i.e., a closed, non-selfintersecting surface, the number of intersections of each ray with the model is even. Voxels along the ray between an odd and the next even intersection are considered to be inside the model and are thus filled. This process is repeated along the Y and Z axis to minimize errors. For now, the internal density distribution is assumed to be uniform, i.e., all voxels will assume the same value. Future work will include non-uniform density distributions. Second, the scanner geometry and properties have to be known to facilitate an exact comparison between measured and simulated radiographs. To simulate a cone-beam geometry the ASTRAtoolbox requires the source and detector position as well as the detector and pixel size, as shown in Fig. 7. The sample object is assumed to be centered in the coordinate system origin. The majority of industrial X-ray radiograph scanners are line scanners, i.e. the detector size is 1 pixel in the direction of movement. During simulation the sample remains in the origin but source and detector move in the opposite conveyor direction over a distance
M. van Dael et al. / Postharvest Biology and Technology 128 (2017) 33–43
37
Fig. 6. The pear shape model fitted to 3 different point clouds, with the color scale indicating the model deformation from the mean pear shape in mm. Below, bar charts display the weights of the first 15 components expressed in standard deviations. These 15 components explain 99% of the shape model variance.
Fig. 7. Workflow to simulate a radiograph from a shape model. The fitted shape model is converted to a voxel-representation using a ray-intersection algorithm. This is then used to generate a radiograph using the appropriate geometric properties.
equal to the detector pixel size. The resulting 1-pixel wide projections images are stitched together to form a complete radiograph. In this study, the source and detector were assumed to be 500 mm above and 50 mm below the conveyor respectively, while the detector was assumed to have a width of 200 mm. To study the effect of pixel size on detection accuracy, radiographs were simulated for pixel sizes of 0.5, 1 and 2 mm. To obtain a radiograph using these parameters, a ray-tracing algorithm is used within the ASTRA toolbox, where a ray is defined as one of the lines connecting the source with an individual detector pixel. The ASTRA toolbox uses a GPU implementation of the slice-interpolated projection scheme (Xu and Mueller, 2006; Palenstijn et al., 2011) which uses bilinear interpolation within each slice and integrates the results in a trapezoidal fashion along each ray. 2.3.3. Step 3. Comparison of radiographs Comparing measured and simulated radiographs can be done in several ways. Standardized methods such as the structural similarity index (SSIM, Wang et al., 2004, 2003), histogram comparison (Cao and Petzold, 2006) and Fourier ring correlation (Banterle et al., 2013) are available to assess similarity. In this work, a different approach was used for is simplicity. Because the two radiographs are taken from the same product orientation it is straightforward to simply subtract them to obtain a residual image. From this image one or multiple features can be determined to be used in a classification process. Here, a Bayesian classifier is used
for its simplicity. For every feature, this classifier determines the normal distribution in both classes from the training set. New samples are assigned to the class where it has the largest a posteriori score, i.e. the chance that the feature has the given value considering the given distribution. When multiple features are considered, naïve Bayesian classifiers assign the sample to the class with the highest mean a posteriori score. Differences between the radiographs will be present as positive or negative pixel values in the residual image. Some differences might be present here that are not linked to the presence of any defects due to errors in the algorithm, mostly due to misalignment and resampling differences during voxelization. It makes sense to choose a feature that results in a net difference measure so that small, random errors due to the reason described above and sensor noise which are both zerocentered are canceled out. The simplest feature possessing this property is the sum over all residual pixels, and is used in what follows as classification feature. This feature is fed into a 10-fold cross-validated naive Bayesian classifier (Larsen, 2005). This means the training set was randomly subdivided into 10 different groups. 10 different classifiers were trained with 9 of these groups, while each classifier was tested with the omitted group. 2.4. Reference detection methods For both data sets, reference classification methods using only X-ray radiography were used to evaluate the performance of the proposed method.
38
M. van Dael et al. / Postharvest Biology and Technology 128 (2017) 33–43
For the torus dataset, human operators were used as reference since no custom image processing algorithm could be developed that performed above chance level. A GUI was developed that presented the human operators with the 2240 radiographs of torus samples in random order. For each of these the operators could indicate whether they could visually determine if a defect was present or not. This process was repeated for the different detector resolutions. Cavities in pears were segmented from randomly simulated radiographs by applying a local adaptive threshold as shown in Fig. 8. This local adaptive threshold (Davies, 2012) smooths the image by applying a mean filter in a 5-pixel circular sliding window. The resulting image is subtracted from the original and the resulting difference image is Otsu-thresholded. False positives generated by the sample edge were removed with a mask generated by dilating the Otsu-thresholded (Otsu, 1979) projection (Fig. 8c and d) with a 3-pixel radius circular structuring element. After removing pixels without neighbors, the summed area of the remaining blobs was fed into a 10-fold cross-validated naive Bayesian classifier. 3. Results 3.1. Output of the proposed method Fig. 9 shows measured radiographs, and residual images for samples of the torus dataset. In the columns the size of defect varies from 1 mm to 7 mm while in the rows the fractional density of the defect varies from 0% to 80%. Note that while large voids are quite easy to distinguish in most measured radiographs, this becomes increasingly more difficult when the defects become smaller and with a smaller deviation in density, while also being dependent on the position of the defect in the sample and the sample pose relative to the source-detector axis. The main difficulty is distinguishing contrast differences introduced by the general shape of the sample and the presence of potential defects. The strength of the proposed method is that the contrast differences introduced by sample shape that are contained in the simulated radiograph are removed, preserving only those contrast differences introduced by defects. The results are the residual images, revealing defects that are not distinguishable by the naked eye as shown in Fig. 9. The same holds true for the pear dataset as shown in Fig. 10 in which the samples in the two right columns contain defects while those in the left two columns do not. Running the complete algorithm, including model fitting and classification takes on average less than a second per sample. The larger and hollow defects make them relatively easy to see in the measured radiographs as shown in column a. Column b shows the output of the reference detection method as described in Section 2.4,
while the residuals resulting from the proposed method are shown in column c. 3.2. Classification accuracy 3.2.1. Torus shapes Confusion matrices for classification of the torus dataset using the proposed and reference methods for detector pixel sizes of 0.5, 1 and 2 mm are shown in Table 1. These confusion matrices are computed for the entire dataset, i.e., for the full range of defect sizes and densities. This shows that the proposed method results in better overall accuracy in comparison to the reference method, i.e., a single human operator, for all detector resolutions. This is also shown in Fig. 11, which shows the average classification accuracy in function of defect size for detector resolutions of 0.5, 1 and 2 mm. The proposed method results in a slightly increased false positives rate (1.3–1.6%) while the reference method results in a larger number of false negatives (23.7–50.3%). In practice this means that the proposed method will reject a small number of samples without defects, while using the reference method a larger number of samples with defects will remain undetected. Comparing the proposed to the reference method, this results in a 0.6, 3.4 and 4% loss of negative predictive value (NPV) with increasing detector pixel size. However, there is a considerably larger gain in positive predictive value (PPV) of 13.1, 7.5 and 11.4% respectively for the different detector pixel sizes. Thus with the proposed method less samples that have been classified as ‘without defects’ will effectively contain them. This is favorable, since the end product will be of higher average quality, e.g. for the 0.5 mm detector pixel size 74.2% of samples classified as without defects using the proposed method will effectively be without defects, while this is only true for 61.1% using the reference method. These values hold for the large range of defect size and density considered. A more detailed view of the power of the proposed method is shown in Fig. 12, displaying the correct classification rate of defects with varying size and density of the defect for both the proposed and reference classification at detector resolutions of 0.5, 1 and 2 mm. At the smallest detector pixel size (0.5 mm), the new method is 100% accurate for defects equal or larger than 3.5 mm radius (0.18 cm3), regardless of their contrast with respect to the object’s material. In case of cavities (defect density 0%), all defects larger than 1.5 mm radius (<0.02 cm3) are correctly detected, for all shapes and orientations of the tori. Classification success becomes smaller with less resolved detectors, but the new method is always better than the reference. The decrease in NPV is caused by a lower average accuracy on the bottom matrix row (samples without defects), while the increase in PPV is caused by an accuracy increase in the rest of the matrix (samples with defects) for the proposed method in comparison to the reference method. Large voids are easier to detect than small defects with a limited
Fig. 8. Steps in the reference defect detection method for pears using X-ray radiography. (a) shows a random radiograph, (b) shows the result of the local adaptive threshold, (c) shows the edge mask, (d) is the result of subtracting (c) from (b) and shows the segmented defect.
M. van Dael et al. / Postharvest Biology and Technology 128 (2017) 33–43
39
Fig. 9. Measured radiographs and resulting residual images for the torus dataset for defect sizes ranging from 1 to 7 mm (columns) and defect densities ranging from 0 to 80%.
Fig. 10. Measured radiographs (a), detected defects using classical methods (b) and resulting residual images for pears with (right two columns) and without (left two columns) defects.
deviation in density, as can be expected intuitively. The reference method accuracy for defects with a density variation of 10% does not rise above 50% for defects of any size, while the proposed method accuracy is 100% for defects larger than 4 mm at 1 mm resolution and 3.5 mm at 0.5 mm resolution. This shows that the proposed method is more sensitive to defects with small density deviations than the reference method. The minimal defect size for a minimal accuracy of 50% is more or less equal for both methods when considering voids. The smaller the density deviation of the defect becomes, the faster the minimal
defect size increases for the reference method. Additionally, the proposed method accuracy for defects larger that this minimal size is stable, while occasional misclassifications still occur with the reference method. This can probably be contributed to some inconsistency introduced by the use of a human operator. 3.2.2. Pears with core breakdown The pears with cavities had a smaller range of defects, both in terms of size as well as density (only voids). The confusion matrices for classification of the pear dataset using the proposed and
40
M. van Dael et al. / Postharvest Biology and Technology 128 (2017) 33–43
Table 1 Confusion matrices showing the classification results of the torus dataset in% as obtained by the proposed method (top) and a human operator as reference (bottom) for X-ray detector sizes of 0.5, 1 and 2 mm. Truth Reference
Proposed No defect
Defect
No defect
Defect
Predicted 0.5 mm
No defect Defect Accuracy PPV NPV
36.2 1.3 86.2 74.2 97.4
12.6 50
37.2 0.3 76.1 61.1 98.2
23.7 38.8
Predicted 1 mm
No defect Defect Accuracy PPV NPV
36.1 1.5 78.6 64.5 96.6
19.9 42.6
37.5 0 71.1 57 100
28.3 34.2
Predicted 2 mm
No defect Defect Accuracy PPV NPV
35.9 1.6 67.8 54 95.2
30.6 31.9
31.4 0.1 49.6 42.6 99.2
50.3 12.2
reference methods are shown in Table 2. It shows that the proposed and reference method reach comparable accuracy on the pear dataset. 4. Discussion Results show that the proposed method performs at least as good as reference methods for large and easy to detect defects such as cavities and considerably better for smaller, more difficult defects with less contrast to the sound tissue. The latter situation is most common in horticultural products. Examples are watercore (Kim and Schatzki, 2000) and Braeburn browing disorder (Herremans et al., 2014) in apples and mealiness in pears (Muziri et al., 2016) where small density variations are caused by locally displaced water and tissue breakdown. The main accuracy gain lies in an increased positive predictive value (PPV), meaning that the end product batch after classification will on average be of higher quality. Improvements of the method are possible and listed below. The pear dataset proves that the method is able to cope with samples that show some limited variation in spatial internal density distribution in the sound tissue, but further work is needed to be able to cope with a larger variation. Examples are apples which show an extensive vascular
structure in CT scans. This could be implemented by using a reference density distribution linked to the shape model parameters or basis functions (Ashburner and Friston, 1999) to assign voxel values in the voxelized reference sample. Misclassification using the proposed method might be related to an inaccurate reference sample, often caused by a problem during shape and pose estimation. The current method can be sensitive to local minima and might have problems dealing with occlusions in the point cloud. Various improved methods have been described in literature though: Schneider and Eisert (2009) combined the pose and shape estimation by linearizing the rigid transformation and solving both together in a least squares sense. Arellano and Dahyot (2012) proposed a mean shift algorithm that uses a Bayesian framework to perform the registration. Future work will include one of these methods, since point cloud data often contain occlusions. The reference method for the torus dataset consisted of a single human operator visually inspecting the measured radiographs for the presence of defects. Only comparing with a single human operator might not be ideal, but in practice every image would be inspected only by a single operator as well. Additionally, for every defect size/density combination at least 10 different images were generated and assessed. Nonetheless it would be interesting in future work to include more human operators to assess operator variability. Both methods falsely show small, low intensity features in the samples without defects but a well-trained classifier omits these during classification. For the reference detection method these can be removed by fine tuning the detection algorithm, a very timeconsuming and application specific process. For the proposed method, these false positives most probably arise from small errors during reference sample estimation and more importantly from local density variations in healthy samples. Since these errors are present in the training set as well, a well-trained classifier makes sure that these samples will still be classified as not containing any defects. The reason the proposed method does not outperform the reference method for the pear dataset can be contributed to the fact that the defects are rather large and completely void, thus easily visible on the radiographs. The results are still valuable since they show that the proposed method can handle samples with a high shape variability without any loss in performance. The proposed method clearly can handle biological samples with a spatially varying internal density without problems. Additionally, it should be noted that when the defects are smaller and show less deviation in density with respect to the rest of the sample, the
Fig. 11. Result of the classifier for the fictive torus dataset for detector resolutions of 0.5, 1 and 2 mm. Full lines show the results of the proposed method, dashed lines show the results of the reference method.
M. van Dael et al. / Postharvest Biology and Technology 128 (2017) 33–43
41
Fig. 12. Correct classification rates on the torus dataset using the proposed method (left) and by a human operator (right) for detector resolutions of 0.5 (top), 1 (middle) and 2 mm (bottom). Colors indicate the classification accuracy for a given defect size – density combination. N = 2240.
Table 2 Confusion matrices showing the classification results in% as obtained with the proposed (left) and reference method (right). Truth Proposed
Predicted
No defect Defect Accuracy PPV NPV
Reference
No defect
Defect
No defect
Defect
70.3 2.1 95.2 96.2 92.2
2.8 24.8
78.3 1.7 95 96 90.3
3.3 16.7
reference method’s performance will start to decline while the proposed method will remain accurate as previously discussed and shown in Fig. 12. The current simulation of X-ray radiographs as implemented in the ASTRA toolbox is purely geometrical, and comparison between simulated and measured radiographs is easy since the whole setup is simulated. To ensure this comparison is still possible in real life, hardware setups the reference sample voxel values will have to be calibrated by for example scanning a reference sample. If this proves inadequate, more advanced X-ray radiograph simulators can be implemented, for example as described in Dhaene et al. (2015). These types of radiograph simulators completely model the physical interactions of the X-rays with physical matter but require detailed information such as X-ray source spectra and detector response curves and are generally far slower than the type of simulators used in this document, resulting in a trade-off between accuracy and speed. The fact that the proposed method detects deviations from an ideal is a big advantage. This implicates that no specific image
processing algorithms need to be developed for different types of defects, as described in Section 2.4. Detection algorithms for various defects and properties have been developed, for example as described in Donis-González et al. (2016, 2015) where very large sets of image features are calculated and fed into advanced classification algorithms. While these methods work, the extensive algorithm development needed for every new application greatly reduces flexibility and increases costs. This highly specialized algorithm development is not needed in the proposed method, greatly increasing flexibility and decreasing costs once a shape model of a sample type has been developed. Although run times (less than a second per sample) are currently still too slow for some practical applications – it is worth noting that the current implementation is completely unoptimized for high-speed performance. Considerable speed improvements are achievable, especially during voxelization. 5. Conclusions A novel method for internal quality inspection based on a combination of 3-D vision and X-ray radiographs is proposed. The proposed method can be implemented on industry-standard X-ray radiograph scanners by adding a 3-D-vision system resulting in much lower costs than for computed tomography (CT) systems while showing significant improvements over classical X-ray radiograph based systems. The method was validated in silico using two datasets and compared to reference methods: a human operator and an algorithm based on classical image processing methods. Results show that the proposed method is more accurate than the reference methods in detecting small defects with a limited density deviation and on par with the reference methods in detecting large cavities. Also, the proposed method has a higher
42
M. van Dael et al. / Postharvest Biology and Technology 128 (2017) 33–43
positive predictive value, meaning the end product is on average of higher quality. Additionally, it should be noted not many improvements can be made to the reference methods, while several improvements have been proposed to increase the proposed method’s performance. Combined with the fact that the proposed method detects deviations from the ideal and is not focused on detecting specific defects or features this shows that the proposed method is a very promising technology. Acknowledgments We would like to thank the Institute for the Promotion of Innovation through Science and Technology (IWT-Vlaanderen, grant numbers 093469 and SBO12033), Flanders Fund for Scientific Research (FWO Vlaanderen, project G.0603.08, project G. A108.10N), the Hercules foundation (project AKUL 09/001), and the EC (project PicknPack FP7-311987) for financial support. The opinions expressed in this document do by no means reflect their official opinion or that of its representatives. References Arellano, C., Dahyot, R., 2012. Shape model fitting algorithm without point correspondence. Proceedings of the 20th European Signal Processing Conference (EUSIPCO), IEEE, pp. 934–938. Ashburner, J., Friston, K.J., 1999. Nonlinear spatial normalization using basis functions. Hum. Brain Mapp. 7, 254–266. doi:http://dx.doi.org/10.1002/(SICI) 1097-0193(1999)7:4<254:AID-HBM4>3.0.CO;2-G. Banterle, N., Bui, K.H., Lemke, E.A., Beck, M., 2013. Fourier ring correlation as a resolution criterion for super-resolution microscopy. J. Struct. Biol. 183, 363– 367. doi:http://dx.doi.org/10.1016/j.jsb.2013.05.004. Barrie Smith, N., Webb, A., 2010. Introduction to Medical Imaging: Physics, Engineering and Clinical Applications. Cambridge University Press. Besl, P.J., McKay, H.D., 1992. A method for registration of 3-D shapes. IEEE Trans. Pattern Anal. Mach. Intell. 14, 239–256. doi:http://dx.doi.org/10.1109/34.121791. Bobelyn, E., Serban, A.-S., Nicu, M., Lammertyn, J., Nicolai, B.M., Saeys, W., 2010. Postharvest quality of apple predicted by NIR-spectroscopy: study of the effect of biological variability on spectra and model performance. Postharvest Biol. Technol. 55, 133–143. doi:http://dx.doi.org/10.1016/j.postharvbio.2009.09.006. Cao, Y., Petzold, L., 2006. Accuracy limitations and the measurement of errors in the stochastic simulation of chemically reacting systems. J. Comput. Phys. 212, 6– 24. doi:http://dx.doi.org/10.1016/j.jcp.2005.06.012. Clark, C.J., Richardson, A.C., Marsh, K.B., 1999. Quantitative magnetic resonance imaging of satsuma mandarin fruit during growth. HortScience 1071–1075. Danckaers, F., Huysmans, T., van Dael, M., Verboven, P., Nicolai, B., Sijbers, J., 2015. Building a statistical shape model of the apple from corresponded surfaces. Chem. Eng. Trans. 44, 49–54. doi:http://dx.doi.org/10.3303/CET1544009. Davies, E.R., 2012. Computer and Machine Vision: Theory, Algorithms, Practicalities. Elsevier. Dhaene, J., Pauwels, E., De Schryver, T., De Muynck, A., Dierick, M., Van Hoorebeke, L., 2015. A realistic projection simulator for laboratory based X-ray micro-CT. Nucl. Instrum. Methods Phys. Res. Sect. B: Beam Interact. Mater. Atoms. 342, 170–178. doi:http://dx.doi.org/10.1016/j.nimb.2014.09.033. Donis-González, I.R., Guyer, D.E., Chen, R., Pease, A., 2015. Evaluation of undesirable fibrous tissue in processing carrots using computed tomography (CT) and structural fiber biochemistry. J. Food Eng. 153, 108–116. doi:http://dx.doi.org/ 10.1016/j.jfoodeng.2014.12.012. Donis-González, I.R., Guyer, D.E., Pease, A., 2016. Postharvest noninvasive assessment of undesirable fibrous tissue in fresh processing carrots using computer tomography images. J. Food Eng. 190, 154–166. doi:http://dx.doi.org/ 10.1016/j.jfoodeng.2016.06.024. Franck, C., Lammertyn, J., Ho, Q.T., Verboven, P., Verlinden, B., Nicolaï, B.M., 2007. Browning disorders in pear fruit. Postharvest Biol. Technol. 43, 1–13. doi:http:// dx.doi.org/10.1016/j.postharvbio.2006.08.008. Haff, R.P., Slaughter, D.C., Sarig, Y., Kader, A., 2006. X-ray assessment of translucency in pineapple. J. Food Process. Preserv. 30, 527–533. doi:http://dx.doi.org/ 10.1111/j.1745-4549.2006.00086.x. Haff, R.P., Saranwong, S., Thanapase, W., Janhiran, A., Kasemsumran, S., Kawano, S., 2013. Automatic image analysis and spot classification for detection of fruit fly infestation in hyperspectral images of mangoes. Postharvest Biol. Technol. 86, 23–28. doi:http://dx.doi.org/10.1016/j.postharvbio.2013.06.003. Hansen, J.D., Schlaman, D.W., Haff, R.P., Yee, W.L., 2005. Potential postharvest use of radiography to detect internal pests in deciduous tree fruits. J. Entomol. Sci. 40 (3), 255–262. Herremans, E., Melado-Herreros, A., Defraeye, T., Verlinden, B., Hertog, M., Verboven, P., Val, J., Fernández-Valle, M.E., Bongaers, E., Estrade, P., Wevers, M., Barreiro, P., Nicolaï, B.M., 2014. Comparison of X-ray CT and MRI of watercore disorder of different apple cultivars. Postharvest Biol. Technol. 87, 42–50. doi: http://dx.doi.org/10.1016/j.postharvbio.2013.08.008.
Kim, S., Schatzki, T.F., 2000. Apple watercore sorting system using X-ray imagery: I. Algorithm development. Trans. ASAE 43, 1695–1702. doi:http://dx.doi.org/ 10.13031/2013.3070. Kotwaliwale, N., Singh, K., Kalne, A., Jha, S.N., Seth, N., Kar, A., 2014. X-ray imaging methods for internal quality evaluation of agricultural produce. J. Food Sci. Technol. 51, 1–15. doi:http://dx.doi.org/10.1007/s13197-011-0485-y. López-Maestresalas, A., Keresztes, J.C., Goodarzi, M., Arazuri, S., Jarén, C., Saeys, W., 2016. Non-destructive detection of blackspot in potatoes by Vis-NIR and SWIR hyperspectral imaging. Food Control 70, 229–241. doi:http://dx.doi.org/ 10.1016/j.foodcont.2016.06.001. Lammertyn, J., Dresselaers, T., Van Hecke, P., Jancsók, P., Wevers, M., Nicolaï, B., 2003a. MRI and x-ray CT study of spatial distribution of core breakdown in conference pears. Magn. Reson. Imaging 21, 805–815. doi:http://dx.doi.org/ 10.1016/S0730-725X(03)00105-X. Lammertyn, J., Dresselaers, T., Van Hecke, P., Jancsók, P., Wevers, M., Nicolaı̈, B., 2003b. Analysis of the time course of core breakdown in conference pears by means of MRI and X-ray CT. Postharvest Biol. Technol. 29, 19–28. doi:http://dx. doi.org/10.1016/S0925-5214(02)00212-0. Larsen, K., 2005. Generalized naive bayes classifiers. ACM SIGKDD Explor. Newsl. 7, 76–81. doi:http://dx.doi.org/10.1145/1089815.1089826. Magwaza, L.S., Opara, U.L., 2014. Investigating non-destructive quantification and characterization of pomegranate fruit internal structure using X-ray computed tomography. Postharvest Biol. Technol. 95, 1–6. doi:http://dx.doi.org/10.1016/j. postharvbio.2014.03.014. Magwaza, L.S., Opara, U.L., Terry, L.A., Landahl, S., Cronje, P.J., Nieuwoudt, H., Mouazen, A.M., Saeys, W., Nicolaï, B.M., 2012. Prediction of Nules Clementine mandarin susceptibility to rind breakdown disorder using Vis/NIR spectroscopy. Postharvest Biol. Technol. 74, 1–10. doi:http://dx.doi.org/10.1016/ j.postharvbio.2012.06.007. Mukundan, R., Ramakrishnan, K., 1998. Moment Functions in Image Analysis— Theory and Applications. World Scientific doi:http://dx.doi.org/10.1142/3838. Muziri, T., Theron, K.I., Cantre, D., Wang, Z., Verboven, P., Nicolai, B.M., Crouch, E.M., 2016. Microstructure analysis and detection of mealiness in Forelle pear (Pyrus communis L.) by means of X-ray computed tomography. Postharvest Biol. Technol. 120, 145–156. doi:http://dx.doi.org/10.1016/j.postharvbio.2016.06.006. Nicolaï, B.M., Beullens, K., Bobelyn, E., Peirs, A., Saeys, W., Theron, K.I., Lammertyn, J., 2007. Nondestructive measurement of fruit and vegetable quality by means of NIR spectroscopy: a review. Postharvest Biol. Technol. 46, 99–118. doi:http://dx. doi.org/10.1016/j.postharvbio.2007.06.024. Nicolaï, B.M., Defraeye, T., De Ketelaere, B., Herremans, E., Hertog, M.L.A.T.M., Saeys, W., Torricelli, A., Vandendriessche, T., Verboven, P., 2014. Nondestructive measurement of fruit and vegetable quality. Annu. Rev. Food Sci. Technol. 5, 285–312. doi:http://dx.doi.org/10.1146/annurev-food-030713-092410. Otsu, N., 1979. A threshold selection method from gray-level histograms. IEEE Trans. Syst. Man Cybern. 9, 62–66. doi:http://dx.doi.org/10.1109/TSMC.1979.4310076. Palenstijn, W.J., Batenburg, K.J., Sijbers, J., 2011. Performance improvements for iterative electron tomography reconstruction using graphics processing units (GPUs). J. Struct. Biol. 176, 250–253. doi:http://dx.doi.org/10.1016/j. jsb.2011.07.017. Palenstijn, W.J., Batenburg, K.J., Sijbers, J., 2013. The ASTRA tomography toolbox. 13th International Conference on Computational and Mathematical Methods in Science and Engineering, CMMSE 2013 1139–1145. Patil, S., Ravi, B., 2005. Voxel-based representation display and thickness analysis of intricate shapes. Ninth International Conference on Computer Aided Design and Computer Graphics (CAD-CG’05), IEEE, pp. 415–422. doi:http://dx.doi.org/ 10.1109/CAD-CG.2005.86. Rogge, S., Defraeye, T., Herremans, E., Verboven, P., Nicolaï, B.M., 2015. A 3D contour based geometrical model generator for complex-shaped horticultural products. J. Food Eng. 157, 24–32. doi:http://dx.doi.org/10.1016/j.jfoodeng.2015.02.006. Schatzki, T.F., Haff, R.P., Young, R., Can, I., Le, L.C., Toyofuku, N., 1996. Defect detection in apples by means of X-ray imaging. In: Meyer, G.E., DeShazer, J.A. (Eds.), Photonics East ’96 International Society for Optics and Photonics, , pp. 176–185. doi:http://dx.doi.org/10.1117/12.262857. Schneider, D.C., Eisert, P., 2009. Fitting a morphable model to pose and shape of a point cloud. Vision, Modeling, and Visualization Workshop 2009. Proceedings 93–100. Shewfelt, R., Brueckner, B., Florkowski, S.E.P.J., 2012. Postharvest Handling: A Systems Approach. Academic Press. van Aarle, W., Palenstijn, W.J., De Beenhouwer, J., Altantzis, T., Bals, S., Batenburg, K. J., Sijbers, J., 2015. The ASTRA Toolbox: a platform for advanced algorithm development in electron tomography. Ultramicroscopy 157, 35–47. doi:http:// dx.doi.org/10.1016/j.ultramic.2015.05.002. van Dael, M., Lebotsa, S., Herremans, E., Verboven, P., Sijbers, J., Opara, U.L., Cronje, P. J., Nicolaï, B.M., 2016a. A segmentation and classification algorithm for online detection of internal disorders in citrus using X-ray radiographs. Postharvest Biol. Technol. 112, 205–214. doi:http://dx.doi.org/10.1016/j. postharvbio.2015.09.020. van Dael, M., Verboven, P., Dhaene, J., Van Hoorebeke, L., Sijbers, J., Nicolaı̈, B., 2016. Automated quality control and selection. PCT/EP2016/055718. Verboven, P., Kerckhofs, G., Mebatsion, H.K., Ho, Q.T., Temst, K., Wevers, M., Cloetens, P., Nicolaï, B.M., 2008. Three-dimensional gas exchange pathways in pome fruit characterized by synchrotron x-ray computed tomography. Plant Physiol. 147, 518–527. doi:http://dx.doi.org/10.1104/pp.108.118935. Wang, Z., Simoncelli, E.P., Bovik, A.C., 2003. Multiscale structural similarity for image quality assessment. The Thrity-Seventh Asilomar Conference on Signals,
M. van Dael et al. / Postharvest Biology and Technology 128 (2017) 33–43 Systems & Computers, IEEE, pp. 1398–1402. doi:http://dx.doi.org/10.1109/ acssc.2003.1292216. Wang, Z., Bovik, A.C., Sheikh, H.R., Simoncelli, E.P., 2004. Image quality assessment: from error visibility to structural similarity. IEEE Trans. Image Process. 13, 600– 612. doi:http://dx.doi.org/10.1109/TIP.2003.819861. Xu, Fang, Mueller, K., 2006. A comparative study of popular interpolation and integration methods for use in computed tomography. 3rd IEEE International
43
Symposium on Biomedical Imaging: Macro to Nano, 2006, IEEE, pp. 125–1255. doi:http://dx.doi.org/10.1109/ISBI.2006.1625152. Zhang, L., McCarthy, M.J., 2013. Assessment of pomegranate postharvest quality using nuclear magnetic resonance. Postharvest Biol. Technol. 77, 59–66. doi: http://dx.doi.org/10.1016/j.postharvbio.2012.11.006.