- Email: [email protected]
Automatic crack monitoring using photogrammetry and image processing
Measurement 46 (2013) 433–441
Contents lists available at SciVerse ScienceDirect
Measurement journal homepage: www.elsevier.com/locate/measurement
Automatic crack monitoring using photogrammetry and image processing J. Valença a,d,⇑, D. Dias-da-Costa b,e, E. Júlio a,f, H. Araújo c, H. Costa a,d a
ICIST, Av. Rovisto Pais, 1049-001 Lisboa, Portugal INESC, Coimbra, Rua Antero de Quental 199, 3000-033 Coimbra, Portugal ISR, Department of Electrical and Computer Engineering, University of Coimbra, Rua Luís Reis Santos, 3030-290 Coimbra, Portugal d Department of Civil Engineering, Institute Polytechnic of Coimbra, Rua Pedro Nunes – Quinta da Nora, 3030-199 Coimbra, Portugal e Department of Civil Engineering, University of Coimbra, Rua Luís Reis Santos, 3030-788 Coimbra, Portugal f Department of Civil Engineering, Instituto Superior Técnico, Technical University of Lisbon, Av. Rovisco Pais, 1049-001 Lisbon, Portugal b c
a r t i c l e
i n f o
Article history: Received 28 July 2011 Received in revised form 27 June 2012 Accepted 27 July 2012 Available online 7 August 2012 Keywords: Laboratorial tests Crack monitoring Crack characterisation Image processing Photogrammetry Monitoring
a b s t r a c t This manuscript presents an integrated approach for automatic crack monitoring combining photogrammetry and image processing. In summary, the strain field obtained from photogrammetric data is used to map the cracked areas where image processing is applied. All processing is completely automatic since only a threshold value, related to the width of the crack, needs to be provided. Direct Shear Tests (DSTs) have been selected for calibration, validation and also as an experimental example. In conclusion, critical areas, the corresponding crack pattern and all related measures (e.g. crack width, length, area or path) could be provided for any stage of loading, until the complete failure of the specimens. Furthermore, all outputs require low computational cost, thus allowing monitoring vast campaigns of laboratorial tests. Ó 2012 Elsevier Ltd. All rights reserved.
1. Introduction The development of automatic procedures for monitoring the most relevant parameters during laboratorial tests is of fundamental importance. In the case of concrete members, the evolution of the crack pattern can assist a complete understanding of the overall structural behaviour up to failure. Most existing approaches to obtain this pattern are still hand-sketched based and the crack opening is often evaluated using measuring magnifiers or crack width rulers. Different authors have recently applied image processing techniques aiming to characterise cracks [1–6]. These algorithms are based on a sudden variation of intensity of pixels to detect a crack. If other sources of discontinuities are present on the surface of the specimen (e.g. voids, ⇑ Corresponding author at: Department of Civil Engineering, Institute Polytechnic of Coimbra, Rua Pedro Nunes – Quinta da Nora, 3030-199 Coimbra, Portugal. E-mail address: [email protected] (J. Valença). 0263-2241/$ - see front matter Ó 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.measurement.2012.07.019
stains, shadows), a sudden variation of intensity of pixels is also identified by the algorithm and the detection of a real crack becomes compromised. Therefore, all existing approaches are still applied to small areas under very controlled conditions. A new approach has been proposed in [7] to automatically monitor the crack pattern for a complete surface. However, the selection of critical regions of analysis still depends on the user, which precludes an automatic procedure. In order to circumvent these drawbacks a new method is herein presented. It was specifically developed to automatically monitor crack evolution, on a given surface, by combining photogrammetry and image processing to minimise user intervention. The strain field, directly computed from photogrammetric data, is used to map the cracked areas where image processing is applied. Furthermore, the spatial resolution is also computed thus allowing measuring the width, length or area of any selected crack. This procedure has been calibrated and validated using Direct Shear Tests (DSTs). One of these specimens is herein analysed in detail.
J. Valença et al. / Measurement 46 (2013) 433–441
The manuscript is organised as follows: the experimental set-up is briefly described in Section 2. In Section 3, the most relevant steps from the proposed approach are provided and thoroughly detailed. In Section 4, the method is illustrated with a DST test. Finally, in Section 5, the most important conclusions are summarised. 2. Experimental set-up 2.1. DST specimens
1200 1000
Load (kN)
434
800 600 400 200
The tested DST specimens were made up of three layers, each with 100 250 200 mm3. The middle layer was cast 28 days after the others, with a 50 mm gap at the bottom (see Fig. 1). The interfaces between layers were crossed by four 6 mm steel connectors. The concrete compressive strength, fcm,cube, reached on average 50 MPa at 28 days of age and the Young’s modulus, Ecm, reached 34 GPa. According to Eurocode-2 [8], this leads to an expected average tensile strength, fctm, of 3.2 MPa. The specimens were loaded by enforcing a vertical displacement of the lower plate of the testing machine (see Fig. 1). The typical load vs. average vertical displacement curve is shown in Fig. 2. The red circular dots represent the three stages selected to illustrate the method, described in detail in Section 3. 2.2. Image acquisition Homogeneous and diffused light pattern conditions were created at the laboratory to prevent discontinuities from being introduced in the images by varying light conditions. A low-cost digital camera, Nikon D200, was installed on a tripod at approximately 90 cm from the surface of the specimen. A remote shutter allowed to maintain the camera stabilised during the complete test. All images were acquired at maximum resolution (3872 2592 pixels) using a 28 mm lens. 3. Method for crack monitoring The new method for crack monitoring is composed by the following main steps, herein detailed: (i) surface preparation (Section 3.1); (ii) data processing (Section 3.2); (iii)
0 0
5
10
15
20
Exp. average vertical displacement (mm) Fig. 2. Load vs. displacement curve.
photogrammetry (Section 3.3); and (iv) image processing (Section 3.4). In summary, the proposed procedure for crack monitoring requires prior preparation of the surface of the specimen, consisting of marking a regular grid made up of circular targets according to the precision required (Section 3.1). Then, the whole process is fully automatic and includes the following main sub-steps: computation of the spatial resolution inside the grid and of the corresponding coordinates of each target, at any selected stage of the test, using a single image per stage (Section 3.3.1); estimation of the strain field using finite element procedures to map the critical areas (Section 3.3.2); application of image processing to the critical region and elimination of surface imperfections (which are not consistent with the geometry of a crack) by means of mathematical morphology techniques (Section 3.4). 3.1. Surface preparation The surface of the specimen is first painted with a regular grid of red circular targets. The grid spacing is established according to the precision required for the strain field (see [9]). In this case a 10 10 mm2 grid was selected. 3.2. Data processing Images are acquired using the set-up described in Section 2.2. All images are pre-processed in order to select: (i) the region of interest (ROI) (compare Figs. 1 and 3a); and (ii) the most relevant of the three RGB1 bands (red–R, green–G and blue–B) from the visible spectrum (see Fig. 3b–d). In this case, the ROI is the surface of the specimen, corresponding to 1847 1866 pixels. The colour of the circular targets was defined to allow their easy removal at image processing stage. In fact, if the R-band is selected, the targets almost vanish and cracking can be straightforwardly characterised, whereas the B or G-bands are most suitable for target detection. Additionally, by subtracting the R-band to the G-band, the crack pattern is removed from
Fig. 1. DST specimen.
1 For interpretation of colour in Fig. 3, the reader is referred to the web version of this article.
J. Valença et al. / Measurement 46 (2013) 433–441
435
Fig. 3. Image of the specimen: (a) RGB; (b) R-band; (c) G-band; and (d) B-band.
the image, providing the most adequate way of target detection. Target detection is performed using the Hough transform to identify the geometrical centre of all targets, at all defined stages [10]. In this case the transform is computed by taking the gradient of the processed band and the edge points that lie along the outline of a circle (with a given radius determined by the painted grid) contribute to the transform at the centre of the circle. Therefore, the peaks in the Hough transform provide directly the centres of the circular targets. Each time a peak is detected, points within one-half of the provided radius are excluded as possible new centres to avoid detecting the same target repeatedly. Fig. 5. Average error map in target detection (in mm).
3.3. Photogrammetry This step includes obtaining the plane coordinates of each target at any selected stage of the test (Section 3.3.1). Then, the strain field can be computed and used to define critical areas (Section 3.3.2). 3.3.1. Homography Fig. 4 shows a perspective image of a target on the world plane being projected on the image plane. If all target points of interest lie on the same world plane, a simplified 3D world coordinate system can be chosen such that Z = 0. Consequently, the following correspondence between the centre of the target point, at each image, and the real plane coordinates can be established [11]:
X ¼ xHx;
ð1Þ T
where X = (X, Y, T) are the world plane homogeneous coordinates, x = (x, y, 1)T are the homogeneous coordinates in
Fig. 4. Perspective image of a point on the world plane being projected on the image plane.
Fig. 6. Half-interval at 95% confidence map in target detection (in mm).
the image, H is the 3 3 homography matrix and x is the scale factor. Notice that: (i) the ‘1’ in x means that all projected targets are placed at a finite distance from the camera (more details concerning an homogeneous coordinate system are given in [12]) and (ii) the ‘T’ in X is introduced to represent a point placed on the world plane and on the same projective ray (see Fig. 4) being mapped into the same image point (see also [12]). Therefore, it can be assumed without any loss of generality that T = 1 and Eq. (1) can be explicitly written as:
436
J. Valença et al. / Measurement 46 (2013) 433–441
Fig. 9. Map of the first principal strain at stage 1.
Fig. 7. RMS error map in target detection (in mm).
Fig. 8. Error map for the homography (in mm).
8 9 2 h1 > = 6 Y ¼ x4 h 4 > : > ; 1 h7
h2 h5 h8
38 9 h3 > = 7 h6 5 y : > : > ; 1 h9
ð2Þ
The total number of unknowns of the problem is eight since only the ratio of the matrix elements is significant. Furthermore, the following two equations can be written, per target point from Eq. (2):
x¼
8 < X ¼ h1 xþh2 yþh3 () Xðh7 x þ h8 y þ h9 Þ ¼ h1 x þ h2 y þ h3 1 h7 xþh8 yþh9 ) h7 x þ h8 y þ h9 : Y ¼ h4 xþh5 yþh6 () Yðh7 x þ h8 y þ h9 Þ ¼ h4 x þ h5 y þ h6 h7 xþh8 yþh9
ð3Þ Eq. (3) is applied to each of the ‘n’ available target points leading to the following global system of equations: 0
x1 B B0 B B x2 B B 0 Ah ¼ B B B .. B . B B @ xn 0
y1 0 y2 0 .. . yn 0
1 0 0 x1 1 0 0 x2 .. .. . . 1 0 0 xn
0 y1 0 y2 .. . 0 yn
0 1 0 1 .. . 0 1
x1 X 1 x1 Y 1 x2 X 2 x2 Y 2 .. . xn X n xn Y n
y1 X 1 y1 Y 1 y2 X 2 y2 Y 2 .. . yn X n yn Y n
0 1 1 h1 X 1 B h2 C CB C C Y 1 CB h C CB B 3C X 2 C h C CB CB 4 C C Y 2 CB h5 C ¼ 0: CB C .. CB B h6 C C . CB C B C h7 C C X n AB B C @ h8 A Y n h9 ð4Þ
The latter system of equations can be solved exactly if four targets are used. However, typically there are many more targets and the problem is over determined. Consequently, the system is solved by minimising the norm |Ah|. In this case, it can be shown that the solution is directly given by the eigenvector corresponding to the least eigenvalue of ATA [11,12]. The parameters of the homography H are computed immediately before starting the experimental test. It is highlighted that: (i) the world plane target points are all placed on the surface and (ii) parallelism is not required between world and image planes (see Fig. 4). The real coordinates are given by the grid spacing, considering the horizontal and vertical axes aligned with grid. For the subsequent stages of the test, the computed parameters are used to obtain the real coordinates of all targets at the new positions, thus scaling and orientating all images. Ten images taken immediately before the experimental test have been used to define the error of the homography. This number was considered sufficient taking into account the characteristics of the digital camera and that these pictures were obtained within a short period of time, i.e. under the same lighting conditions. Therefore, these images allow estimating the error for: (i) target detection algorithm (Section 3.2) and (ii) homography (Eq. (1)). The average difference of the coordinates at each target is represented in Fig. 5, whereas the half-interval at 95% confidence is represented in Fig. 6. The upper bound at confidence 95%, evaluated by adding Figs. 5 and 6, led to an average value of 0.034 mm and a maximum value of 0.105 mm. Fig. 7 shows the Root Mean Square (RMS) precision derived from the differences of the coordinates at each target. This error was on average of 0.027 mm, having a minimum of 0.011 mm and a maximum of 0.084 mm. Fig. 8 contains the error map for the homography. This error was on average of 0.241 mm, whereas the minimum and maximum values were, respectively, of 0.022 mm and 0.875 mm. Notice that this is in fact a systematic error, affecting the real coordinates of each pixel. Nevertheless, since only relative displacements are needed, this error tends to vanish. Therefore, in the following sections, the error in the measurements of the displacements is
J. Valença et al. / Measurement 46 (2013) 433–441
Fig. 10. Map of the first principal strain at stage 2.
Fig. 11. Map of the first principal strain at stage 3.
437
3.3.2. Critical areas With the displacement defined at each target using the homography, the corresponding strain field can be evaluated. For that purpose, an auxiliary triangular mesh using the targets is defined by means of a Delaunay triangulation [13]. This ensures that no target is contained in any triangle’s circumscribed circle. Then, the strain-nodal displacement matrix at each target is used to compute the strain field, for all stages of monitoring (see [9] for a more detailed description of the method). The precision of the strain field can be roughly estimated by assuming a uniform local strain field. For the adopted grid, this corresponds to 0.3%, meaning that the method is adequate for detecting large deformations such as those related to the crack localisation process. The maps for the first principal strain are represented in Figs. 9–11, for all stages. The critical areas are the ones that have the first principal strain above a threshold defined by the user. This value must take into account: (i) the strain corresponding to the tensile strength of the material (in this case 0.01%); (ii) the average size of the cracks that the user aims to identify; (iii) the minimum value below which surface discontinuities are incorrectly identified as cracks; and (iv) the average pixel size, below which measurements have no physical meaning. In the case herein considered it was decided to evaluate cracks above 0.1 mm, therefore within the photogrammetric precision reached with the adopted grid. Fig. 12a represents the auxiliary mesh for stage 3 which envelops first principal strains above 1% (see Fig. 11). Since the pixel size was 0.16 mm in average, it is not possible to detect cracks below this value by image processing. Therefore, image processing was constrained to an auxiliary mesh corresponding to a threshold value of 2.5% (see Fig. 12b). 3.4. Image processing
considered to be generated by the target detection algorithm. This value is assumed to be the average of the upper bound at 95% confidence, i.e., 0.034 mm.
Image processing is now constrained to the mesh, composed by the triangles identified in the previous section (for each stage).
Fig. 12. Auxiliary mesh for stage 3 defined with different thresholds: (a) et = 1.0% and (b) et = 2.5%.
438
J. Valença et al. / Measurement 46 (2013) 433–441
Fig. 13. Image processing: (a) definition of a rectangular mask; (b) high-pass filter; (c) binarization; (d) noise reduction; (e) joining operation; (f) filling gaps; (g) shape detection; and (h) final discontinuity inside a triangular element.
Fig. 14. LROI: (a) initial stage; (b) failure stage with correction; and (c) failure stage without correction.
Each triangular element is sequentially processed and a rectangular mask is defined using the maximum and minimum coordinates of the triangle (Fig. 13a). This rectangular mask is required for applying image processing algorithms. At the end, all data outside the triangle are disregarded. A high-pass filter is applied to enhance the image (Fig. 13b) and Otsu’s method [14] is applied to binarize the image (Fig. 13c). A set of mathematical morphology [15,16] operations is then applied to the resultant rectangular image aiming to: (i) reduce the noise due to the binarization process by removing areas representing less than
10 pixels (Fig. 13d); (ii) join pixels within a neighbourhood of three pixels (Fig. 13e); (iii) filling gaps (Fig. 13f); and (iv) removing round objects (Fig. 13g). The last operation allows eliminating surface discontinuities (e.g. voids), which due to having circular geometry, cannot be cracks. This is performed by computing the parameter 4p(Area/Perimeter2) and removing all shapes above 0.75. Other discontinuities (cracks or imperfections) which do not met this criterion, such as linear discontinuities, are not removed from the analysis. The discontinuities contained in the triangular element are stored (Fig. 13h). The complete crack pattern results
439
J. Valença et al. / Measurement 46 (2013) 433–441
Fig. 15. Crack width rulers placed in front of the specimen.
from adding up all stored data. In order to characterise any crack, the user has to select a local region of interest (LROI). Then, the following sequence of three steps allows evaluating the length, width and area: (1) localisation of the crack on the surface; (2) definition of the boundary edges; and (3) evaluation of geometrical properties of the crack (see [7] for a more detailed explanation). In order to ensure that the selected LROI remains the same along the several stages of the analysis, its position changes according to the displacement of the boundary, evaluated using the target points. Neglecting this correction can lead to errors (see Fig. 14a to c). Crack width rulers (CWR) were placed in front of the specimen immediately before starting the test (see Fig. 15). Each ruler is composed by a set of lines with increasing thickness, ranging from 0.1 to 4.0 mm. The CWR were printed in a plotter capable of defini minimum thickness lines of 0.042 mm and precision of ±0.2%. All CWR have been measured using the described procedure. Fig. 16a contains the measured widths vs. known values showing a coefficient of correlation near 1.0. It should be highlighted that, due to the resolution of the images, the thinnest lines (0.1 and 0.2 mm) cannot always be measured since these are within the pixel size. For measurements above 0.2 mm, the accuracy is always below 0.1 mm (see Fig. 16b).
Fig. 17. Crack pattern at stage 1.
4. Results In this section, the main results obtained in relation to crack pattern and crack characterisation are summarised. 4.1. Crack pattern Figs. 17–19 present a crack pattern which is above 0.25 mm wide on average, obtained from image processing. The coloured lines represented for each target are scaled in order to represent the average crack opening expected from photogrammetry. Three scales have been chosen: 0.1–0.25 mm, 0.25–0.5 mm, and above 0.5 mm. It is also remarked that these lines are perpendicular to the first principal strain, which can be directly correlated with the direction of the crack opening. Stage 1 was taken immediately after the first sudden drop of the applied load (see Fig. 2), which corresponds to the first crack appearing on the surface of the specimen. This crack can be seen in Fig. 17. Between stages 1 and 2,
0.10
y = 0.9964x R² = 0.996
0.08
3.0
accuracy (mm)
measured width (mm)
4.0
2.0
0.06 0.04
1.0 0.02 0.0 0.0
1.0
2.0
3.0
reference width (mm)
4.0
0.00 0.0
1.0
2.0
3.0
reference width (mm)
Fig. 16. Measured widths vs. known values of all CWRs: (a) correlation and (b) accuracy.
4.0
440
J. Valença et al. / Measurement 46 (2013) 433–441
Fig. 18. Crack pattern at stage 2.
Fig. 20. Otsu’s method applied to the complete surface of the specimen.
Time (in min.) 03:01
04:04
05:22
Fig. 21. Evolution of the crack profile.
Area / Length
8
Fig. 19. Crack pattern at stage 3.
1.00
area (mm 2 ) length (mm) width (mm)
0.75
6 0.50 4 0.25
2
the opening of the left interface takes place. The sliding occurring at the interface can be identified by the direction of the first principal strain. Simultaneously, there is a propagation of the crack placed at the right layer of the specimen towards the right interface (Fig. 18). In the last stage, there is an increase of the number of cracks present at the bottom of the right layer which is related to crushing (Fig. 19). It should also be mentioned that the crack detected at the left interface still exhibits noise due to the presence of imperfections on the surface. However, with the proposed approach, in which image processing is focused on critical areas, results improve significantly. This can be observed by comparing the crack pattern obtained with the direct application of the Otsu’s method, after band-R
0 00:00
Width
10
01:06
02:12
03:18
04:24
0.00 05:30
Time (min) Fig. 22. Geometrical properties of the crack contained by the selected lroi.
enhancement, to the three layers of the specimen at stage 3 (see Figs. 19 and 20). 4.2. Crack characterisation Fig. 21 shows the evolution of the bottom of the right crack (shaded area) in the corresponding selected LROI.
J. Valença et al. / Measurement 46 (2013) 433–441
Fig. 22 contains a graphical representation of the evolution of the geometrical properties of the discontinuity along several stages. Note that: (i) the crack was only detected four minutes after the start of the test: (ii) the average crack opening is similar to the one obtained from photogrammetry (see Figs. 17–19). 5. Conclusions To adequately assess the evolution of the crack pattern of concrete members is quite valuable for the complete understanding of the behaviour of concrete structures. The method herein proposed, is fully capable of monitoring a whole concrete surface from crack onset up to failure. This is performed by using a combined approach based on photogrammetry and image processing. Moreover, all processing is completely automatic since only a threshold value, related to the width of the crack, needs to be provided. All outputs require low computational cost and are fast to obtain, allowing the method to be used for monitoring vast campaigns of laboratorial tests. It should also be mentioned that the surface of the specimen does not require a specific treatment, since only a regular grid of circular targets needs to be painted. Note that this grid can in fact be used additionally to obtain more valuable information as, for instance, both displacement and strain fields (see [9,17]). The major drawback of the proposed approach is related to the initial discontinuities, e.g. imperfections on the concrete surface, that can be incorrectly detected as cracks along the test. This was specially relevant in the specimen studied, since it was made up of three concrete layers, cast at different ages, resulting in significant initial discontinuities at the concrete joints. The authors are currently addressing this issue, by studying and developing a different technique aiming at overcoming this difficulty. Acknowledgements This work is supported by FEDER funds through the Operational Programme for Competitiveness Factors –
441
COMPETE – and by Portuguese funds through FCT – Portuguese Foundation for Science and Technology under Project No. FCOMP-01-0124-FEDER-020275 (FCT ref. PTDC/ ECM/119214/2010). References [1] M. Sachtleber, Z. Zhao, D. Raabe, Experimental investigation of plastic grain interaction, Materials Science and Engineering A 336 (2002) 81–87. [2] H. Thomas, S. Cantré, Applications of low-budget photogrammetry in the geotechnical laboratory, The Photogrammetric Record 24 (2009) 332–350. [3] P. Dare, H. Hanley, C. Fraser, B. Ridel, W. Niemeier, An operational application of the automatic feature extraction: the measurement of cracks in concrete structures, The Photogrammetric Record 17 (2002) 453–464. [4] S. Sinha, P. Fieguth, Segmentation of buried concrete pipe images, Automation in Construction 15 (2006) 47–57. [5] T. Yamaguchi, S. Nakamura, R. Saegusa, S. Hashimoto, Image-based crack detection for real concrete surfaces, IEEJ Transactions on Electrical and Electronic Engineering 3 (2008) 128–135. [6] L. Barazzetti, M. Scaioni, Crack measurement: development, testing and applications of an automatic image-based algorithm, ISPRS Journal of Photogrammetry and Remote Sensing 64 (2009) 285–296. [7] J. Valença, D. Dias-da-Costa, E.N.B.S. Júlio, Characterisation of concrete cracking during laboratorial tests using image processing, Construction and Building Materials 28 (2012) 607–615. [8] CEN, EN 1992-1-1, Eurocode 2: Design of Concrete Structures - Part 1-1: General Rules and Rules for Buildings, European Committee for Standardization (CEN), 2004. [9] D. Dias-da-Costa, J. Valença, E. Júlio, Laboratorial test monitoring applying photogrammetric post-processing procedures to surface displacements, Measurement 44 (2011) 527–538. [10] D. Ballard, Generalizing the hough transform to find arbitrary shapes, Pattern Recognition 13 (1981) 111–122. [11] A. Criminisi, I. Reid, A. Zisserman, Single view metrology, International Journal Computers Vision 40 (2000) 123–148. [12] R. Hartley, A. Zisserman, Multiple View Geometry in Computer Vision, second ed., Cambridge University Press, 2003. ISBN: 0521540518. [13] C.B. Barber, D.P. Dobkin, H. Huhdanpaa, The quickhull algorithm for convex hulls, ACM Transactions on Mathematical Software 22 (1996) 469–483. [14] N. Otsu, A threshold selection method from gray-level histogram, IEEE Transactions on System Man Cybernetics, SMC-9 (1979) 62–66. [15] R. Gonzales, R. Woods, Digital Image Processing, second ed., Prentice Hall, 2002. [16] S. Marchand-Maillet, Y. Sharaiha, Binary Digital Image Processing, A Discret Approach, Academic press, London, 2000. [17] J. Valença, E. Júlio, H. Araújo, Application of photogrammetry to structural assessment, Experimental Techniques, Wiley 36 (2012).
Contents lists available at SciVerse ScienceDirect
Measurement journal homepage: www.elsevier.com/locate/measurement
Automatic crack monitoring using photogrammetry and image processing J. Valença a,d,⇑, D. Dias-da-Costa b,e, E. Júlio a,f, H. Araújo c, H. Costa a,d a
ICIST, Av. Rovisto Pais, 1049-001 Lisboa, Portugal INESC, Coimbra, Rua Antero de Quental 199, 3000-033 Coimbra, Portugal ISR, Department of Electrical and Computer Engineering, University of Coimbra, Rua Luís Reis Santos, 3030-290 Coimbra, Portugal d Department of Civil Engineering, Institute Polytechnic of Coimbra, Rua Pedro Nunes – Quinta da Nora, 3030-199 Coimbra, Portugal e Department of Civil Engineering, University of Coimbra, Rua Luís Reis Santos, 3030-788 Coimbra, Portugal f Department of Civil Engineering, Instituto Superior Técnico, Technical University of Lisbon, Av. Rovisco Pais, 1049-001 Lisbon, Portugal b c
a r t i c l e
i n f o
Article history: Received 28 July 2011 Received in revised form 27 June 2012 Accepted 27 July 2012 Available online 7 August 2012 Keywords: Laboratorial tests Crack monitoring Crack characterisation Image processing Photogrammetry Monitoring
a b s t r a c t This manuscript presents an integrated approach for automatic crack monitoring combining photogrammetry and image processing. In summary, the strain field obtained from photogrammetric data is used to map the cracked areas where image processing is applied. All processing is completely automatic since only a threshold value, related to the width of the crack, needs to be provided. Direct Shear Tests (DSTs) have been selected for calibration, validation and also as an experimental example. In conclusion, critical areas, the corresponding crack pattern and all related measures (e.g. crack width, length, area or path) could be provided for any stage of loading, until the complete failure of the specimens. Furthermore, all outputs require low computational cost, thus allowing monitoring vast campaigns of laboratorial tests. Ó 2012 Elsevier Ltd. All rights reserved.
1. Introduction The development of automatic procedures for monitoring the most relevant parameters during laboratorial tests is of fundamental importance. In the case of concrete members, the evolution of the crack pattern can assist a complete understanding of the overall structural behaviour up to failure. Most existing approaches to obtain this pattern are still hand-sketched based and the crack opening is often evaluated using measuring magnifiers or crack width rulers. Different authors have recently applied image processing techniques aiming to characterise cracks [1–6]. These algorithms are based on a sudden variation of intensity of pixels to detect a crack. If other sources of discontinuities are present on the surface of the specimen (e.g. voids, ⇑ Corresponding author at: Department of Civil Engineering, Institute Polytechnic of Coimbra, Rua Pedro Nunes – Quinta da Nora, 3030-199 Coimbra, Portugal. E-mail address: [email protected] (J. Valença). 0263-2241/$ - see front matter Ó 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.measurement.2012.07.019
stains, shadows), a sudden variation of intensity of pixels is also identified by the algorithm and the detection of a real crack becomes compromised. Therefore, all existing approaches are still applied to small areas under very controlled conditions. A new approach has been proposed in [7] to automatically monitor the crack pattern for a complete surface. However, the selection of critical regions of analysis still depends on the user, which precludes an automatic procedure. In order to circumvent these drawbacks a new method is herein presented. It was specifically developed to automatically monitor crack evolution, on a given surface, by combining photogrammetry and image processing to minimise user intervention. The strain field, directly computed from photogrammetric data, is used to map the cracked areas where image processing is applied. Furthermore, the spatial resolution is also computed thus allowing measuring the width, length or area of any selected crack. This procedure has been calibrated and validated using Direct Shear Tests (DSTs). One of these specimens is herein analysed in detail.
J. Valença et al. / Measurement 46 (2013) 433–441
The manuscript is organised as follows: the experimental set-up is briefly described in Section 2. In Section 3, the most relevant steps from the proposed approach are provided and thoroughly detailed. In Section 4, the method is illustrated with a DST test. Finally, in Section 5, the most important conclusions are summarised. 2. Experimental set-up 2.1. DST specimens
1200 1000
Load (kN)
434
800 600 400 200
The tested DST specimens were made up of three layers, each with 100 250 200 mm3. The middle layer was cast 28 days after the others, with a 50 mm gap at the bottom (see Fig. 1). The interfaces between layers were crossed by four 6 mm steel connectors. The concrete compressive strength, fcm,cube, reached on average 50 MPa at 28 days of age and the Young’s modulus, Ecm, reached 34 GPa. According to Eurocode-2 [8], this leads to an expected average tensile strength, fctm, of 3.2 MPa. The specimens were loaded by enforcing a vertical displacement of the lower plate of the testing machine (see Fig. 1). The typical load vs. average vertical displacement curve is shown in Fig. 2. The red circular dots represent the three stages selected to illustrate the method, described in detail in Section 3. 2.2. Image acquisition Homogeneous and diffused light pattern conditions were created at the laboratory to prevent discontinuities from being introduced in the images by varying light conditions. A low-cost digital camera, Nikon D200, was installed on a tripod at approximately 90 cm from the surface of the specimen. A remote shutter allowed to maintain the camera stabilised during the complete test. All images were acquired at maximum resolution (3872 2592 pixels) using a 28 mm lens. 3. Method for crack monitoring The new method for crack monitoring is composed by the following main steps, herein detailed: (i) surface preparation (Section 3.1); (ii) data processing (Section 3.2); (iii)
0 0
5
10
15
20
Exp. average vertical displacement (mm) Fig. 2. Load vs. displacement curve.
photogrammetry (Section 3.3); and (iv) image processing (Section 3.4). In summary, the proposed procedure for crack monitoring requires prior preparation of the surface of the specimen, consisting of marking a regular grid made up of circular targets according to the precision required (Section 3.1). Then, the whole process is fully automatic and includes the following main sub-steps: computation of the spatial resolution inside the grid and of the corresponding coordinates of each target, at any selected stage of the test, using a single image per stage (Section 3.3.1); estimation of the strain field using finite element procedures to map the critical areas (Section 3.3.2); application of image processing to the critical region and elimination of surface imperfections (which are not consistent with the geometry of a crack) by means of mathematical morphology techniques (Section 3.4). 3.1. Surface preparation The surface of the specimen is first painted with a regular grid of red circular targets. The grid spacing is established according to the precision required for the strain field (see [9]). In this case a 10 10 mm2 grid was selected. 3.2. Data processing Images are acquired using the set-up described in Section 2.2. All images are pre-processed in order to select: (i) the region of interest (ROI) (compare Figs. 1 and 3a); and (ii) the most relevant of the three RGB1 bands (red–R, green–G and blue–B) from the visible spectrum (see Fig. 3b–d). In this case, the ROI is the surface of the specimen, corresponding to 1847 1866 pixels. The colour of the circular targets was defined to allow their easy removal at image processing stage. In fact, if the R-band is selected, the targets almost vanish and cracking can be straightforwardly characterised, whereas the B or G-bands are most suitable for target detection. Additionally, by subtracting the R-band to the G-band, the crack pattern is removed from
Fig. 1. DST specimen.
1 For interpretation of colour in Fig. 3, the reader is referred to the web version of this article.
J. Valença et al. / Measurement 46 (2013) 433–441
435
Fig. 3. Image of the specimen: (a) RGB; (b) R-band; (c) G-band; and (d) B-band.
the image, providing the most adequate way of target detection. Target detection is performed using the Hough transform to identify the geometrical centre of all targets, at all defined stages [10]. In this case the transform is computed by taking the gradient of the processed band and the edge points that lie along the outline of a circle (with a given radius determined by the painted grid) contribute to the transform at the centre of the circle. Therefore, the peaks in the Hough transform provide directly the centres of the circular targets. Each time a peak is detected, points within one-half of the provided radius are excluded as possible new centres to avoid detecting the same target repeatedly. Fig. 5. Average error map in target detection (in mm).
3.3. Photogrammetry This step includes obtaining the plane coordinates of each target at any selected stage of the test (Section 3.3.1). Then, the strain field can be computed and used to define critical areas (Section 3.3.2). 3.3.1. Homography Fig. 4 shows a perspective image of a target on the world plane being projected on the image plane. If all target points of interest lie on the same world plane, a simplified 3D world coordinate system can be chosen such that Z = 0. Consequently, the following correspondence between the centre of the target point, at each image, and the real plane coordinates can be established [11]:
X ¼ xHx;
ð1Þ T
where X = (X, Y, T) are the world plane homogeneous coordinates, x = (x, y, 1)T are the homogeneous coordinates in
Fig. 4. Perspective image of a point on the world plane being projected on the image plane.
Fig. 6. Half-interval at 95% confidence map in target detection (in mm).
the image, H is the 3 3 homography matrix and x is the scale factor. Notice that: (i) the ‘1’ in x means that all projected targets are placed at a finite distance from the camera (more details concerning an homogeneous coordinate system are given in [12]) and (ii) the ‘T’ in X is introduced to represent a point placed on the world plane and on the same projective ray (see Fig. 4) being mapped into the same image point (see also [12]). Therefore, it can be assumed without any loss of generality that T = 1 and Eq. (1) can be explicitly written as:
436
J. Valença et al. / Measurement 46 (2013) 433–441
Fig. 9. Map of the first principal strain at stage 1.
Fig. 7. RMS error map in target detection (in mm).
Fig. 8. Error map for the homography (in mm).
8 9 2 h1 >
h2 h5 h8
38 9 h3 >
ð2Þ
The total number of unknowns of the problem is eight since only the ratio of the matrix elements is significant. Furthermore, the following two equations can be written, per target point from Eq. (2):
x¼
8 < X ¼ h1 xþh2 yþh3 () Xðh7 x þ h8 y þ h9 Þ ¼ h1 x þ h2 y þ h3 1 h7 xþh8 yþh9 ) h7 x þ h8 y þ h9 : Y ¼ h4 xþh5 yþh6 () Yðh7 x þ h8 y þ h9 Þ ¼ h4 x þ h5 y þ h6 h7 xþh8 yþh9
ð3Þ Eq. (3) is applied to each of the ‘n’ available target points leading to the following global system of equations: 0
x1 B B0 B B x2 B B 0 Ah ¼ B B B .. B . B B @ xn 0
y1 0 y2 0 .. . yn 0
1 0 0 x1 1 0 0 x2 .. .. . . 1 0 0 xn
0 y1 0 y2 .. . 0 yn
0 1 0 1 .. . 0 1
x1 X 1 x1 Y 1 x2 X 2 x2 Y 2 .. . xn X n xn Y n
y1 X 1 y1 Y 1 y2 X 2 y2 Y 2 .. . yn X n yn Y n
0 1 1 h1 X 1 B h2 C CB C C Y 1 CB h C CB B 3C X 2 C h C CB CB 4 C C Y 2 CB h5 C ¼ 0: CB C .. CB B h6 C C . CB C B C h7 C C X n AB B C @ h8 A Y n h9 ð4Þ
The latter system of equations can be solved exactly if four targets are used. However, typically there are many more targets and the problem is over determined. Consequently, the system is solved by minimising the norm |Ah|. In this case, it can be shown that the solution is directly given by the eigenvector corresponding to the least eigenvalue of ATA [11,12]. The parameters of the homography H are computed immediately before starting the experimental test. It is highlighted that: (i) the world plane target points are all placed on the surface and (ii) parallelism is not required between world and image planes (see Fig. 4). The real coordinates are given by the grid spacing, considering the horizontal and vertical axes aligned with grid. For the subsequent stages of the test, the computed parameters are used to obtain the real coordinates of all targets at the new positions, thus scaling and orientating all images. Ten images taken immediately before the experimental test have been used to define the error of the homography. This number was considered sufficient taking into account the characteristics of the digital camera and that these pictures were obtained within a short period of time, i.e. under the same lighting conditions. Therefore, these images allow estimating the error for: (i) target detection algorithm (Section 3.2) and (ii) homography (Eq. (1)). The average difference of the coordinates at each target is represented in Fig. 5, whereas the half-interval at 95% confidence is represented in Fig. 6. The upper bound at confidence 95%, evaluated by adding Figs. 5 and 6, led to an average value of 0.034 mm and a maximum value of 0.105 mm. Fig. 7 shows the Root Mean Square (RMS) precision derived from the differences of the coordinates at each target. This error was on average of 0.027 mm, having a minimum of 0.011 mm and a maximum of 0.084 mm. Fig. 8 contains the error map for the homography. This error was on average of 0.241 mm, whereas the minimum and maximum values were, respectively, of 0.022 mm and 0.875 mm. Notice that this is in fact a systematic error, affecting the real coordinates of each pixel. Nevertheless, since only relative displacements are needed, this error tends to vanish. Therefore, in the following sections, the error in the measurements of the displacements is
J. Valença et al. / Measurement 46 (2013) 433–441
Fig. 10. Map of the first principal strain at stage 2.
Fig. 11. Map of the first principal strain at stage 3.
437
3.3.2. Critical areas With the displacement defined at each target using the homography, the corresponding strain field can be evaluated. For that purpose, an auxiliary triangular mesh using the targets is defined by means of a Delaunay triangulation [13]. This ensures that no target is contained in any triangle’s circumscribed circle. Then, the strain-nodal displacement matrix at each target is used to compute the strain field, for all stages of monitoring (see [9] for a more detailed description of the method). The precision of the strain field can be roughly estimated by assuming a uniform local strain field. For the adopted grid, this corresponds to 0.3%, meaning that the method is adequate for detecting large deformations such as those related to the crack localisation process. The maps for the first principal strain are represented in Figs. 9–11, for all stages. The critical areas are the ones that have the first principal strain above a threshold defined by the user. This value must take into account: (i) the strain corresponding to the tensile strength of the material (in this case 0.01%); (ii) the average size of the cracks that the user aims to identify; (iii) the minimum value below which surface discontinuities are incorrectly identified as cracks; and (iv) the average pixel size, below which measurements have no physical meaning. In the case herein considered it was decided to evaluate cracks above 0.1 mm, therefore within the photogrammetric precision reached with the adopted grid. Fig. 12a represents the auxiliary mesh for stage 3 which envelops first principal strains above 1% (see Fig. 11). Since the pixel size was 0.16 mm in average, it is not possible to detect cracks below this value by image processing. Therefore, image processing was constrained to an auxiliary mesh corresponding to a threshold value of 2.5% (see Fig. 12b). 3.4. Image processing
considered to be generated by the target detection algorithm. This value is assumed to be the average of the upper bound at 95% confidence, i.e., 0.034 mm.
Image processing is now constrained to the mesh, composed by the triangles identified in the previous section (for each stage).
Fig. 12. Auxiliary mesh for stage 3 defined with different thresholds: (a) et = 1.0% and (b) et = 2.5%.
438
J. Valença et al. / Measurement 46 (2013) 433–441
Fig. 13. Image processing: (a) definition of a rectangular mask; (b) high-pass filter; (c) binarization; (d) noise reduction; (e) joining operation; (f) filling gaps; (g) shape detection; and (h) final discontinuity inside a triangular element.
Fig. 14. LROI: (a) initial stage; (b) failure stage with correction; and (c) failure stage without correction.
Each triangular element is sequentially processed and a rectangular mask is defined using the maximum and minimum coordinates of the triangle (Fig. 13a). This rectangular mask is required for applying image processing algorithms. At the end, all data outside the triangle are disregarded. A high-pass filter is applied to enhance the image (Fig. 13b) and Otsu’s method [14] is applied to binarize the image (Fig. 13c). A set of mathematical morphology [15,16] operations is then applied to the resultant rectangular image aiming to: (i) reduce the noise due to the binarization process by removing areas representing less than
10 pixels (Fig. 13d); (ii) join pixels within a neighbourhood of three pixels (Fig. 13e); (iii) filling gaps (Fig. 13f); and (iv) removing round objects (Fig. 13g). The last operation allows eliminating surface discontinuities (e.g. voids), which due to having circular geometry, cannot be cracks. This is performed by computing the parameter 4p(Area/Perimeter2) and removing all shapes above 0.75. Other discontinuities (cracks or imperfections) which do not met this criterion, such as linear discontinuities, are not removed from the analysis. The discontinuities contained in the triangular element are stored (Fig. 13h). The complete crack pattern results
439
J. Valença et al. / Measurement 46 (2013) 433–441
Fig. 15. Crack width rulers placed in front of the specimen.
from adding up all stored data. In order to characterise any crack, the user has to select a local region of interest (LROI). Then, the following sequence of three steps allows evaluating the length, width and area: (1) localisation of the crack on the surface; (2) definition of the boundary edges; and (3) evaluation of geometrical properties of the crack (see [7] for a more detailed explanation). In order to ensure that the selected LROI remains the same along the several stages of the analysis, its position changes according to the displacement of the boundary, evaluated using the target points. Neglecting this correction can lead to errors (see Fig. 14a to c). Crack width rulers (CWR) were placed in front of the specimen immediately before starting the test (see Fig. 15). Each ruler is composed by a set of lines with increasing thickness, ranging from 0.1 to 4.0 mm. The CWR were printed in a plotter capable of defini minimum thickness lines of 0.042 mm and precision of ±0.2%. All CWR have been measured using the described procedure. Fig. 16a contains the measured widths vs. known values showing a coefficient of correlation near 1.0. It should be highlighted that, due to the resolution of the images, the thinnest lines (0.1 and 0.2 mm) cannot always be measured since these are within the pixel size. For measurements above 0.2 mm, the accuracy is always below 0.1 mm (see Fig. 16b).
Fig. 17. Crack pattern at stage 1.
4. Results In this section, the main results obtained in relation to crack pattern and crack characterisation are summarised. 4.1. Crack pattern Figs. 17–19 present a crack pattern which is above 0.25 mm wide on average, obtained from image processing. The coloured lines represented for each target are scaled in order to represent the average crack opening expected from photogrammetry. Three scales have been chosen: 0.1–0.25 mm, 0.25–0.5 mm, and above 0.5 mm. It is also remarked that these lines are perpendicular to the first principal strain, which can be directly correlated with the direction of the crack opening. Stage 1 was taken immediately after the first sudden drop of the applied load (see Fig. 2), which corresponds to the first crack appearing on the surface of the specimen. This crack can be seen in Fig. 17. Between stages 1 and 2,
0.10
y = 0.9964x R² = 0.996
0.08
3.0
accuracy (mm)
measured width (mm)
4.0
2.0
0.06 0.04
1.0 0.02 0.0 0.0
1.0
2.0
3.0
reference width (mm)
4.0
0.00 0.0
1.0
2.0
3.0
reference width (mm)
Fig. 16. Measured widths vs. known values of all CWRs: (a) correlation and (b) accuracy.
4.0
440
J. Valença et al. / Measurement 46 (2013) 433–441
Fig. 18. Crack pattern at stage 2.
Fig. 20. Otsu’s method applied to the complete surface of the specimen.
Time (in min.) 03:01
04:04
05:22
Fig. 21. Evolution of the crack profile.
Area / Length
8
Fig. 19. Crack pattern at stage 3.
1.00
area (mm 2 ) length (mm) width (mm)
0.75
6 0.50 4 0.25
2
the opening of the left interface takes place. The sliding occurring at the interface can be identified by the direction of the first principal strain. Simultaneously, there is a propagation of the crack placed at the right layer of the specimen towards the right interface (Fig. 18). In the last stage, there is an increase of the number of cracks present at the bottom of the right layer which is related to crushing (Fig. 19). It should also be mentioned that the crack detected at the left interface still exhibits noise due to the presence of imperfections on the surface. However, with the proposed approach, in which image processing is focused on critical areas, results improve significantly. This can be observed by comparing the crack pattern obtained with the direct application of the Otsu’s method, after band-R
0 00:00
Width
10
01:06
02:12
03:18
04:24
0.00 05:30
Time (min) Fig. 22. Geometrical properties of the crack contained by the selected lroi.
enhancement, to the three layers of the specimen at stage 3 (see Figs. 19 and 20). 4.2. Crack characterisation Fig. 21 shows the evolution of the bottom of the right crack (shaded area) in the corresponding selected LROI.
J. Valença et al. / Measurement 46 (2013) 433–441
Fig. 22 contains a graphical representation of the evolution of the geometrical properties of the discontinuity along several stages. Note that: (i) the crack was only detected four minutes after the start of the test: (ii) the average crack opening is similar to the one obtained from photogrammetry (see Figs. 17–19). 5. Conclusions To adequately assess the evolution of the crack pattern of concrete members is quite valuable for the complete understanding of the behaviour of concrete structures. The method herein proposed, is fully capable of monitoring a whole concrete surface from crack onset up to failure. This is performed by using a combined approach based on photogrammetry and image processing. Moreover, all processing is completely automatic since only a threshold value, related to the width of the crack, needs to be provided. All outputs require low computational cost and are fast to obtain, allowing the method to be used for monitoring vast campaigns of laboratorial tests. It should also be mentioned that the surface of the specimen does not require a specific treatment, since only a regular grid of circular targets needs to be painted. Note that this grid can in fact be used additionally to obtain more valuable information as, for instance, both displacement and strain fields (see [9,17]). The major drawback of the proposed approach is related to the initial discontinuities, e.g. imperfections on the concrete surface, that can be incorrectly detected as cracks along the test. This was specially relevant in the specimen studied, since it was made up of three concrete layers, cast at different ages, resulting in significant initial discontinuities at the concrete joints. The authors are currently addressing this issue, by studying and developing a different technique aiming at overcoming this difficulty. Acknowledgements This work is supported by FEDER funds through the Operational Programme for Competitiveness Factors –
441
COMPETE – and by Portuguese funds through FCT – Portuguese Foundation for Science and Technology under Project No. FCOMP-01-0124-FEDER-020275 (FCT ref. PTDC/ ECM/119214/2010). References [1] M. Sachtleber, Z. Zhao, D. Raabe, Experimental investigation of plastic grain interaction, Materials Science and Engineering A 336 (2002) 81–87. [2] H. Thomas, S. Cantré, Applications of low-budget photogrammetry in the geotechnical laboratory, The Photogrammetric Record 24 (2009) 332–350. [3] P. Dare, H. Hanley, C. Fraser, B. Ridel, W. Niemeier, An operational application of the automatic feature extraction: the measurement of cracks in concrete structures, The Photogrammetric Record 17 (2002) 453–464. [4] S. Sinha, P. Fieguth, Segmentation of buried concrete pipe images, Automation in Construction 15 (2006) 47–57. [5] T. Yamaguchi, S. Nakamura, R. Saegusa, S. Hashimoto, Image-based crack detection for real concrete surfaces, IEEJ Transactions on Electrical and Electronic Engineering 3 (2008) 128–135. [6] L. Barazzetti, M. Scaioni, Crack measurement: development, testing and applications of an automatic image-based algorithm, ISPRS Journal of Photogrammetry and Remote Sensing 64 (2009) 285–296. [7] J. Valença, D. Dias-da-Costa, E.N.B.S. Júlio, Characterisation of concrete cracking during laboratorial tests using image processing, Construction and Building Materials 28 (2012) 607–615. [8] CEN, EN 1992-1-1, Eurocode 2: Design of Concrete Structures - Part 1-1: General Rules and Rules for Buildings, European Committee for Standardization (CEN), 2004. [9] D. Dias-da-Costa, J. Valença, E. Júlio, Laboratorial test monitoring applying photogrammetric post-processing procedures to surface displacements, Measurement 44 (2011) 527–538. [10] D. Ballard, Generalizing the hough transform to find arbitrary shapes, Pattern Recognition 13 (1981) 111–122. [11] A. Criminisi, I. Reid, A. Zisserman, Single view metrology, International Journal Computers Vision 40 (2000) 123–148. [12] R. Hartley, A. Zisserman, Multiple View Geometry in Computer Vision, second ed., Cambridge University Press, 2003. ISBN: 0521540518. [13] C.B. Barber, D.P. Dobkin, H. Huhdanpaa, The quickhull algorithm for convex hulls, ACM Transactions on Mathematical Software 22 (1996) 469–483. [14] N. Otsu, A threshold selection method from gray-level histogram, IEEE Transactions on System Man Cybernetics, SMC-9 (1979) 62–66. [15] R. Gonzales, R. Woods, Digital Image Processing, second ed., Prentice Hall, 2002. [16] S. Marchand-Maillet, Y. Sharaiha, Binary Digital Image Processing, A Discret Approach, Academic press, London, 2000. [17] J. Valença, E. Júlio, H. Araújo, Application of photogrammetry to structural assessment, Experimental Techniques, Wiley 36 (2012).