Application of Fractal Dimension in Industry Practice

Abstract

Today, industrial production lines commonly use off-line and automatic on-line quality monitoring and control. Monitoring and control units analyse data from a production process, and analysis should be able to obtain reliable information that correspond with the character of the data obtained. The character of the data set obtained from production processes or from products can be highly structured in all industrial areas. Structured surface, complex time series (topologically one dimensional signals), difficulty to describe dividing curves are much more common than it can be expected. For this kind of data set, a powerful tool for analysis of complexity —fractal geometry (especially a fractal dimension) should be used. The fractal dimension with a combination of statistical tools is an interesting and powerful tool for complex data quantification, for tracing the source of poor quality, production optimization and investigating the source of instability of production process subsystems in industrial applications. The methodology for evaluation of complex and irregular data was developed and applied in industrial practice. This methodology searches appropriated parameters for a complex evaluation of data. Only the chosen parameters are used for a complete analysis of the data in order to reduce processing time.

1. Introduction

A demand for objective measurement and control methods for materials, processes and production processes stems from continuously increasing pressure from competitors to improve the quality of products. However, description of many complex and irregular structures (e.g. defects, surfaces, cracks, signals from dynamic processes) is almost impossible by conventional methods. The application of fractal geometry, which is successfully used in science, appears to be a powerful approach. The industrial application of the fractal dimension (FD) is generally experimental [1, 2], but application to production is possible and brings benefits. Currently, there are tools to monitor three basic data format types: topological one-dimensional signals (analysis of time series for control systems, roughness from a surface roughness tester), dividing curves (application for the corrugation test of flat glass and surface roughness—iron aluminides in comparison with carbide-nickel steel in contact with molten glass), and digitalized photos (evaluation of 2D pictures of surface defects).

It is possible to compare statistical tools and FD, but the comparison should be done carefully. To claim that the FD is better than statistics and vice versa is impossible, because the FD provides added information about the character of data set description and does not substitute standard tools. Moreover, the FD should not be used separately for the reason that the dimension does not give all the information about data set captures. FD using together with added parameters (statistics, topology, spectral analysis, etc.) brings benefits and is recommended.

2. Description of one-dimensional signals

Dynamic subsystems can be found in many production processes, and they have a strong influence on the production. Data measured from production sensors contain mentioned dynamic effect. Therefore, the time series are structured and the typical statistical data evaluations are not often sufficient.

Also roughness from a surface roughness tester is in many cases complex, and using standard roughness parameters may be not satisfactory.

To analyse one-dimensional signals, we use statistical methods, power spectral analysis, and an estimate of the FD [3, 4]. The estimation of the FD is calculated using the rescaled range method [5] and the box counting method from an “iso-set” [6], also compass counting and EEE method can be used.

2.1. Rescaled range analysis

The rescaled range analysis (R/S) represents method for estimating the FD of self-affine fractals and uses statistical tools. The method is based on an analysis of a changed interval of time series. Consider an interval of time series, length w. Within this interval, one can define two quantities:

R(w), the range taken by the values of y in the interval (vector y—the vertical axis of the time series—represent the time series without a sampling interval x—the horizontal axis of the time series). The range is measured with respect to a trend in the interval, where the trend is estimated simply as difference between maximum and minimum of the interval. This subtracts the average trend in the interval.

S(w), the standard deviation of the first differences Δy of the values of y within the interval. The first differences of the y's are defined as the differences between the values of y at some location x and y at the previous location on the x axis:

$\text{Δ}y\left(x\right)=y\left(x\right)-y\left(x-\text{Δ}x\right)$E1

where Δx is the sampling interval, that is, the interval between two consecutive values of x.

The rescaled range R/S(w_RS) is defined as follows:

$R/S\left(w_{RS}\right)=⟨\frac{R\left(w_{RS}\right)}{S\left(w_{RS}\right)}⟩$E2

where w is the interval length and the angled brackets <R(w_RS)> represent the average of a R(w_RS) value numbers. The self-affinity incurs that one expects the range taken by the values of y in an interval of length w to be proportional to the interval length to a power equal to the Hurst exponent H, that is,

$R/S\left(w_{RS}\right)={\text{w}}_{\text{RS}}^{\text{H}}$E3

For a given window length w_RS, the input series in a number of intervals of length w_RS are then subdivided, R(w_RS) and S(w_RS) in each interval are measured and R/S(w_RS) is calculated as the average ratio R(w_RS)/S(w_RS), see Eq. (1). Mentioned process is then reiterated for several window lengths. Logarithms of R/S(w_RS) are subsequently plotted versus the logarithms of w_RS. Considering the self-affine trace, this plot follows a straight line with a slope equals to the Hurst exponent H. FD of the time series can be calculated from the relationship between the Hurst exponent H and the FD:

$D_{RS}=2-H$E4

where R/S dimension D_RS denotes the FD estimated from the Rescaled Range analysis. The R/S dimension has value from 1 to 2, and Hurst exponent has value from 0 to 1. More about the method can be found in Ref. [5].

2.2. Box dimension from “Iso-set”

An “iso-set” is constructed from the time series as shown in Figure 1 that contains zeros and ones, and the FD is estimated from this set. The “iso-set” (as time series) can be generated by two basic ways, where the ones represent crossings of a pre-selected threshold. The “iso-set” is generated by setting of suitable thresholds and marking the time at which the time series cross these thresholds (Figure 1). The threshold values can be perceptually dependent on the time series average (a floated threshold value, in Figure 1) or can be pre-selected fixed values.

Figure 1.

Construction of iso-set by marking each time at which trace crossed the chosen threshold values and a principle of box dimension computing.

The FD of the “iso-set” is estimated by using box counting method that is described in Section 4.4. The principle of the box dimension method used for “iso-set” is given in Figure 1. (The box size r_B=t_s.b.) Starting from box size r_B = t_s (t_s is sampling time interval), the number of boxes that contain a crossing is recorded. The box size is then increased by the factor b and the procedure continues until the entire “iso-set” is contained in one single box. This is illustrated for the factor b = 2. The box dimension D_B is determined from the central slope of the regression line of the Richardson-Mandelbrot plot (logarithmic dependence between log₂N(r) and log₂r_B). For more information about the method, please see Refs. [6–8].

2.3. Compass counting

The estimated Compass Dimension expresses the degree of complexity of the profile. A compass method [5, 9, 10] is based on measuring of the profile (curve) using different ruler sizes (Figure 2A) according to the equation:

Figure 2.

Estimation of the fractal dimension by the compass method.

$L_i\left(r_i\right)=N_i\left(r_i\right).r_i$E5

where L_i is the length in i-step of the measurement, r_i is the ruler size and N_i is the number of steps needed for the measurement.

If the profile is fractal, and hence the estimated FD is larger than the topological dimension, then the length measured increases as the ruler size is reduced. The logarithmic dependence between log₂N(r_i) and log₂r_i (Richardson-Mandelbrot plot) is shown in Figure 2B. The Compass Dimension is then determined from the slope s of the regression line:

$D_C=1-s=1-\frac{\text{Δ}{\log }_{2}L\left(r\right)}{\text{Δ}{\log }_{2}r}$E6

2.4. Relative and proportional length

The rate of signal (profile) deformation might be evaluated by its relative length L_R. This fast and reliable method measures the ratio of the profile length l_PIXEL (solid line in Figure 2A) using the smallest ruler (1 pixel) r_PIXEL and the length of the projection l (Figure 2A):

$L_R=\frac{l_{PIXEL}}{l}$E7

Another similar approach is to compute the proportional length of the profile L_P. The proportional length is the ratio of the profile length measured with a defined ruler l_r (e.g. dashed line in Figure 2A) and the length measured with the maximum ruler l_{r max} (the length between the first and the last point of the profile):

$L_P=\frac{l_r}{l_{rmax}}$E8

2.5. EEE method

The method is based on length evaluation of a curve (signal). The curve is defined by measured values and they are isolated points x₁, x₂, …,x_n in the range y(x₁), y(x₂), …, y(x_n). The dots represent local extremes (minima and maxima). Unnecessary extremes are classified with a defined rule on the curve, and a new simplified function is defined by the remaining points. For the next classification, the new function is used.

An example of a function defined by points and connected into the linear by parts function f is in Figure 3. The function f is in the first step purged of points which are not local extremes using the rule:

Figure 3.

Function f and its local extremes.

First, the difference proportion of the dependent variables y to the independent variables x between neighbouring points x_i, x_i+1 is determined from function f:

$\text{Δ}f\left(x_i,x_{i+1}\right)=\frac{\text{Δ}y\left(x_i\right)}{\text{Δ}{x}_i}=\frac{y\left(x_{i+1}\right)-y\left(x_i\right)}{x_{i+1}-x_i}$E9

Second, irrelevant points, concretely points where the difference Δf has the same sign, are eliminated. Remained points (marked by black dots in Figure 3) are considered as the local extremes. Extremes are supplemented by first and last points. The prepared points form a simplified function g, Figure 4. A relative length of the function g is computed, and the result is saved. Based on the absolute length of the function the relative length, L_{R 1} (Section 2.4) is evaluated from point to point and divided by the length of its x axis projection.

Figure 4.

Simplified function g and its local extremes.

The elimination of insignificant extremes procedure is used to the simplified function g. The procedure applies functions formed from function g minima and maxima, Figure 5, the functions g_min and g_max are extended by the g function’s first and last points, as it was mentioned above. The function g_max is formed by the g function’s maxima (dotted line in Figure 5). There are local maxima presented in this function by black dots in Figure 5. The definition of local maxima using above-mentioned rule, but only maxima are used (obtained minima are forgotten). The function g_min is generated from function’s minima (dashed line). Obtained local minima are presented as black dots using above-mentioned rule; however, only minima are used (obtained maxima are forgotten).

Figure 5.

Functions g_max and g_min generated from local extremes of function g.

The local minima and maxima of the functions g_max and g_min are used for the generation of the function g_red, Figure 6. In this function, again local minima and maxima are defined using the rule (Figure 6, black dots). These final local extremes of the function g (Figure 4, black dots) and the first and last points from the function g define the function h, Figure 7. The relative length L_{R 2} of function h is computed, and the result is saved.

Figure 6.

Function g_red generated from maxima of function g_max and from minima of function g_min.

Figure 7.

Simplified function h and its local extremes.

The similar procedure is subsequently used for the new simplified function h. Global extremes and the first and last point (Figure 7, black dots) define the function k. The function is formed from the global minimum and maximum of all functions (f, g, h), and therefore the analysis is stopped. All functions are represented in Figure 8.

Figure 8.

Function f and simplified functions g, h, k.

The analysis steps j are plotted versus g, h, k, function’s computed relative lengths L_{R j} (Section 2.4.) Figure 9. The dependence between the relative lengths L_{R j} and steps of elimination j is computed by a sufficient regression function. Concretely, it can be represented by a regression line (Figure 9), a quadratic function or a hyperbolical function. Parameters of those regression functions are used for the function f determination. (Logarithmic axis use is not beneficial.)

Figure 9.

Plot of the relation between steps of analysis and relative length, regression line.

Signal (time series) obtained from simulation of fractional Brownian motion using Cholesky-Factorization of the related covariance matrix (FBM) was used to test the developed method. An example of testing signal is generated using the input Hurst coefficient H = 0.4, in Figure 10. The coefficient represents the character of signal and can achieve value 0 to 1 (lower coefficients generate more complex functions, and further information can be found in Refs. [5, 10]). The dependence between relative lengths and elimination steps is shown in Figure 11.

Figure 10.

Simulation of the time series using fractional Brownian motion.

Figure 11.

Relation between number of steps of analysis and relative length for the time series FBM, H=0.4.

A standard process to estimate the FD using the dependence between log₂L_{R j} and log₂j was tested. However, such representation was not beneficial as in this method the rulers are not used with different lengths. For that reason, the length of rulers cannot be used in the plot.

The regression function, which fits the most to describe a relation of relative lengths L_{R j} and the number of steps j, is estimated by the hyperbolic regression model:

$L_{Rj}=\frac{d}{j+a}+b$E10

This might be calculated by parameters d and a. Parameter b is always set as b = 1. Parameter a needs to be computed numerically with application of an error function.

For verification of the EEE method, over 900 simulated time series from FBM were used with the Hurst coefficient between 0.1 and 0.95. The dependence between Hurst coefficients and the average value of parameters a and d is shown in Figure 12.

Figure 12.

Hurst coefficient versus average value of parameters a and d of hyperbolic function.

EEE dimension for self-affine fractals can be estimated on base of parameter d:

$D_{EEE}=1+\left|d\right|$E11

2.6. Example of time series results

There are numerous possible applications of a fast and accurate description of time series (signals) from production process sensors using fractal geometry and statistical analyses.

Figure 13 represents results of three time series analysis from a glass tank: temperature of a tank main arch. The time series are from one position of production processes in different times. Box and R/S dimensions of the time series are written under the pictures. The complexity is estimated by one number, larger number represents higher complexity. First time series has relatively the largest complexity, and third is relatively the smoother (Figure 13).

Figure 13.

Time series of glass tank main arch temperatures from a first sensor T1 and these box and R/S dimension in different time.

Figure 14 shows results of a tank siege analysis, time series were obtained during defined time, and whole signals were analyzed using the standard deviation, the R/S dimension, and the box dimension. The average temperature indicates a temperature profile of the tank siege and an implicit temperature profile of the molten glass in contact to the tank siege. The box dimensions were computed for threshold values k = {0.2, 0.5,0.7}, see Figure 14. From positions 5 and 6, the FD decreases and between positions 6 and 8 rises. The large dimension indicates a smoother time series, and the higher dimension represents a complex time series. Glass melt is in permanent movement in the tank. This movement is important for a good quality of glass melt in the end of the glass tank. However, this movement is impossible to monitor by standard measurements. The character (complexity) of time series represents changes in temperature. These changes are caused by the movement of molten glass. The decrease of the complexity occurs where a change of longitudinal glass currents is expected. This shows that fractal analysis can be used for the detection of molten glass currents [3, 8] (Figure 14).

Figure 14.

Results of glass tank siege time series analysis.

2.7. Examples of roughness analysis results

The raw data reading (unfiltered reading) from a surface roughness tester presented by a curve is labelled as a profile. Obtained parameters can be divided into three groups:

• parameters of frequency describing surface profile spacing parameter and corrugation frequency characterization (e.g. Sm—Mean Spacing),

• parameters of amplitude describing depth characterization (Ra—Average Roughness, Std—Standard Deviation, Rz—Mean Roughness, Rt—Maximum Roughness, Depth, etc.).

• parameters of complexity and deformation describing FD estimation by Compass Dimension (Dc), by Profile Proportional Length (L_P), EEE method, or Profile Relative Length (L_R).

Mentioned statistic parameters, amplitude, and frequency are widely used in industry. Surface profile parameters, such as Maximum Roughness, Average Roughness, Mean Spacing, and Mean Roughness Depth, are defined by ISO 4287-1997 [12] standard. Complexity and deformation parameters were chosen on the basis of previous experiences.

For better result comparison, the dimension is multiplied by 1000 (D_C1000).

We analysed 14 surfaces produced by five different processes and in different conditions, Table 1. Figure 15 shows 28 samples (with 14 surfaces) [13]. The analysed structures were chosen so as to be different and to cover the most common surfaces in industrial practice. Chosen samples were made purposely from identical material. This allows us to subsequently ignore material properties and to analyse the change of technological parameters and the influence of used technology.

Sample	Technology of surfaces production
1	Polished surface to maximum gloss
2	Ballotini (glass beads) blasting, grain size F120 (mean diameter 0.109 mm)
3	Corundum blasting, grain size F36 (mean diameter 0.525 mm)
4	Corundum blasting, grain size F12 (mean diameter 1.765 mm)
5	Electro-erosion machining 29A
6	Electro-erosion machining 42A
7	Electro-erosion machining 54A
8	Sandpaper, K400
9	Emery cloth, 120
10	Emery cloth, 80
11	Vertical milling machine, milling cutter 20 mm, 120 rpm, feed 30 mm/min
12	Grinding wheel, 98A 60 J 9 V C40
13	Grinding wheel, 96A 36P 5V
14	Vertical milling machine, milling cutter 20 mm, 120 rpm, feed 240 mm/min

Table 1.

Correlation coefficients of selected parameters.

Figure 15.

Analysed samples with machined surfaces.

Measurement was realized on a surface roughness tester Mitutoyo SV 2000. Parameters were set as follows: stylus measuring speed: 0.5 mm/s; positioning: 2 mm/s; traverse range: 50 mm; linearity of traverse: 0.3 µm/50 mm. Standard type of stylus with a 60° angle with a measuring force: 0.75 mN was used. The length of measurement is 4800 µm, and the sampling interval is 0.5 µm.

All samples (two samples with the same surface) were measured in nine positions, each position in three directions, x, y, and transversely. All data obtained in the form of unfiltered profiles were used for analyses. Matlab platform was used for a data evaluation, and necessary software tools were developed.

Based on a linear correlation from obtained parameters (Section 3.1) and simplifying results, we can specify suitable parameters for evaluation of these types of data: Average Roughness, Pa (parameter of amplitude, Figure 16), Mean Spacing, Sm (parameter of frequency Figure 17), Compass Dimension, D_{C 1000} (parameter of complexity and deformation, Figure 18). Diverse information from data can be provided by these three parameters. EEE dimension (parameter of complexity and deformation) is depicted in Figure 19. More information can be found in Ref. [13].

Figure 16.

Results of Pa parameter.

Figure 17.

Results of Sm parameter.

Figure 18.

Results of FD estimation, D_{C 1000}.

Figure 19.

Results of FD estimation, D_{EEE 1000}.

3. Classification of topological one-dimensional dividing curve

The research focuses on the application of the methodology for a quantification of metal surface changes and on an objectification of corrugation test for flat glass. Tools for analysis (the estimation of FD and statistical tools) are in the principle the same as for the topological one-dimensional signals, described in Section 2. The following text shows examples of FD application in practice.

3.1. Example: surface roughness changes after corrosion tests

Quantification of metal surface changes after exposition on air or in a glass melt is important for objective comparison of materials corrosion resistance. Example is showed in relatively new materials: iron aluminides [14], compared with currently used chrome-nickel steels in contact with molten glass. The sample roughness changed during an interaction of metal surface with molten glass, and the effect of disruption can be evaluated after the end of the corrosion test. The obtained roughness of metal surface was quantified using FD and statistical tools.

First, a digital camera (in a light optical microscope) takes a photograph of a metal surface profile from metallographic sample (five photographs from one sample), Figure 20A.

Figure 20.

Image analysis of dividing curves between alloy and glass: A—gray scale image from light optical microscopy, B—evaluated dividing curve between the alloy and glass (boundary curve), C—parameter of amplitude (Rt—maximum roughness) and parameter of frequency (Sm—mean spacing), D—computing of compass dimension, E—compass dimension D_C computed from slope.

Second, a dividing curve is generated from obtained images by a software tool. Software exactly defines a curve between material alloys and its surrounding, Figure 20B. The generated dividing curve is evaluated by the FD (a compass dimension multiplied by 1000, D_{C 1000}), Figure 20D and E. The average standard deviation of all the curves (STD), average mean spacing(Sm), Figure 20C, and the average maximum roughness of all the curves (Rt), Figure 20C, are then described using statistics [15–17].

Examples of analysis are shown on materials: Fe-14Al-6Cr (iron aluminides) and EN X8CrNi25-21 (chrome-nickel steels). More information about tests and materials can be found in Ref. [17]. The dividing curves are depicted in Figure 21.

Figure 21.

Examples of dividing curves chrome-nickel steel material and iron aluminide, after static glass melt effects in different temperatures and results of analyses (compass dimension multiplied by 1000, D_{c 1000}).

The evaluation of roughness parameters was carried out on 10 places for each sample (each alloy and time interval). The corrosion attack of the tested alloy may also be described by the roughness of the surface. Many parameter types can be used for a quantification of the metal roughness. Parameters can be divided into three groups: parameters of amplitude, parameters of frequency, and parameters of complexity and deformation (described in Section 2.7). In this field of research, a filtered profile is not being used. For this reason, the Average Roughness is called Pa, maximum Roughness is denoted Pt etc.

Graphs in Figures 22–24 show the results of analysis for the dividing curves between alloys and glass. Average values were used in order to compare results.

Figure 22.

The average value of maximum roughness of profile, Pt as a function of time for corrosion in molten soda-lime glass.

Figure 23.

The average value of mean spacing of profile, Sm as a function of time for corrosion in molten soda-lime glass.

Figure 24.

The average value of fractal dimension estimation, compass dimension of profile D_{C 1000} as a function of time for corrosion in molten soda-lime glass.

The analysis using the developed methodology has two steps. The first step is a specification of appropriate parameters for fast and reliable analysis for data evaluation. Mentioned methodology contains 22 parameters. However, only chosen parameters were used for a complete analysis in order to simplify the analysis and to reduce the processing time.

Some parameters linearly correlate with the others (they provide similar information about the data). To evaluate the parameters objectively, Pearson's correlation coefficients were computed, see Table 2. A correlation between chosen parameters is clearly visible (Pa and Pt). The L_R parameter correlates less with the parameters Pa and Pt, but still significantly. The D_{C 1000} parameter correlates less with the parameters Pa and Pt, but correlates with the L_R parameter. The Sm parameter does not correlate with other parameters. On the base of a linear correlation from obtained parameters and simplifying results, we can specify suitable parameters for evaluation of these types of data: Maximum Roughness, Pt (parameter of amplitude), Mean Spacing, Sm (parameter of frequency), Compass Dimension, D_{C 1000} (parameter of complexity and deformation). Diverse information from data can be provided by these three parameters

	Relative length, LR[-]	Compass dimension, DC 1000 [-]	Mean spacing, Sm [um]	Average roughness, Pa [um]	Maximum roughness, Pt [um]	Standard deviation, Std [um]
Standard deviation, Std [um]	0.90	0.83	0.50	0.99	0.98	1
Maximum roughness, Pt [um]	0.88	0.85	0.47	0.96	1
Average roughness, Pa [um]	0.91	0.82	0.49	1
Mean spacing, Sm [um]	0.33	0.30	1
Compass dimension, DC 1000 [-]	0.92	1
Relative length, LR[-]	1

Table 2.

List of analysed samples with their production properties.

Second, it is possible to objectively describe a character of metal structure after corrosion attack. We are able to draw the conclusions:

Parameter of amplitude, Maximum Roughness Pt (Figure 22), shows deepness of corrosion attack of metal surface by glass melt. The average dividing curve deepness of alloy Fe-14Al-6Cr grows slowly up to 168 h because iron aluminide dissolved slowly and uniformly in the molten glass than EN X8CrNi25-21 [17–19]. The average dividing curve deepness of austenitic steel grows from raw state during attack to 48 h (Figure 21). After 48 h, it is apparent that corrosion protrusions penetrate less deeply into the surface of steel. It does not mean any increase of corrosion resistance, but probably a progress of corrosion attack is more uniform in this time period. This should be analysed by other methods, for example, the measurement of weight loss, chemical analyses, etc. The maximal deepness of corrosion attack is after 168 h.

Parameter of frequency, Mean Spacing Sm (Figure 23), shows surface profile spacing. The parameter describes corrugation frequency of the dividing curve after corrosion attack, how many wavelets can be observed on the surface. For both alloys grows the parameter up to 72 h, where the maximum is. The frequency is connected to the corrosion mechanism.

Parameter of complexity and deformation, Compass Dimension D_{C 1000} (Figure 24), shows level of dividing curve complexity. In case of alloy Fe-14Al-6Cr, the average complexity of surface expressed as Compass Dimension grows from 24 to 48 h and then slowly decreases. On the other hand, the average surface complexity of alloy EN X8CrNi25-21 increases significantly from raw state during attack up to 48 h, then falls down, and grows to the maximum after 72 h. After 96 h, slight decrease of D_{C 1000} was observed. The significant decrease in value of D_{C 1000} after 48 h is probably related to more uniform a progress of corrosion. For understanding of the steel corrosion mechanism, it is necessary to carry out further analysis.

It is obvious that dividing curves of alloy Fe-14Al-6Cr after corrosion attack are smoother and less complex. It seems that the corrosion resistance of Fe-14Al-6Cr is higher than steel EN X8CrNi25-21. Austenitic steel showed corrosion protrusions due to probably intergranular corrosion (preferential attack of some phases at grain boundaries). However, this statement should be supported by structural and phase analysis.

3.2. Example: corrugation test

The optical test using a zebra plate is widely used measurements for mass production. The test is one of the many important measurements, and it is used in a wide range of situations: as a production control by manufacturers of float glass; as a quality control by glass processors on the glass they buy; as a production control of products (laminated glass, thermal treated glass, etc.); as a quality control by the final customer on the glass they buy.

The corrugation test (Figure 25) is focused on the reflection while another type of tests specialize on the passage of light though the glass (test of distortion). The test is based on the reflection of light off a glass sample sheet from a skew striped plate. The zebra plate is 1 × 2 m with 25 mm wide black strips at an angle of 45°. Glass sheet is laid on a Table 4 m distanced from zebra plate. An observer is distanced additional 4 m from the table. The observer subjectively evaluated the corrugation of reflected light based on a comparison with etalons. The quality of the sheet is classified using a rate dimensionless number from 1.5 to 3.5. The evaluation is carried out off-line in a dark-room. Samples of flat glass are obtained from an on-line production process, and they are cut from the whole width of production glass ribbon. Figure 26 shows good and poor quality of glass sheets during the corrugation test. Using a small angle of observation caused the relatively extreme “distortion”.

Figure 25.

Scheme of corrugation test, A—zebra-plate, B—observer, C—table with window-glass sheet.

Figure 26.

Good and poor quality of glass sheets with a specification of measured parts.

A measuring system, which consists of both hardware and software, has been developed (Figure 27) according to executed experiments. The system is an important intermediate stage for development of the on-line measurement and its calibration. It is also a suitable solution for processors and customers, without an expensive on-line measurement. The hardware of the system includes the zebra-pate (Figure 27A) and a table for the glass sheet sample (C) as in the standard subjective test. The operator is replaced by a scanning unit (digital camera with power adaptor, B), a control unit (PC, D), a connection between the digital camera and the control unit (E), a ball head (F), a camera support system (placed on a top-wall, G), and a system for gripping the table (H) to keep it in a defined place.

Figure 27.

Scheme of measuring system for objective monitoring.

An operator lays a tested glass sheet on the table into the defined place, fills in a form in the Corrugation software and by pressing a button he starts the evaluation. The software initiates the communication with the camera and the image is captured. The image analysis starts after downloading the image from the camera to the computer. In a while, the operator can read the evaluated quality on the PC screen.

The analysis of a scanned image is used for the detection of a glass sheet position, for the necessary rotation and shifting, and for the generation of curves from detected boundaries between the light and dark areas of the sheet in the image. The curves are generated using a thresholding and a detection of contours. They are converted into sequences of points, where the axis X is defined from an ideal reflection, Figure 28. For curve obtaining, two thresholds are used, as two different types of defects on the surface of the bottom part glass sheet can be found: a primary and a secondary corrugation Figure 29. (Different character of the reflected image can be found on heat treated glass and also the top and bottom parts can differ.) The methodology enables to describe both corrugations separately, and it evaluates the impact of the technological production parts to the quality of the production. The obtained curves are then evaluated by software tools: statistic; measuring the curve length; and estimating its FD.

Figure 28.

Conversion of the curve to a sequence of points and its evaluation using range R and relative length L_R.

Figure 29.

Primary and secondary corrugations of bottom part glass sheet.

Our research proved that only one parameter does not enable to describe all types of corrugation. We had to collect limited number of adequate parameters that evaluate completely the quality. An extensive analysis was performed to achieve the best conformity between a subjective evaluations and the evaluated quality using computed parameters. We used hundreds of measurements from one experienced operator with his specification of quality in the current scale. The results showed that three parameters (of the separated curves from the contours) are important to reliably assess the quality of flat glass:

• complexity (smoothness, used the compass dimension mentioned in Section 2.3);

• range of waviness from an ideal line;

• rate of deformation (using the relative length of the curves, Section 2.4).

All parameters are measured in pixels, and the parameters for one sample are defined as average values of all curves. Six quality parameters (the average values from two thresholds) are obtained and used for the final evaluation of the quality. The six quality parameters for each bottom part of a glass sheet obtained from the corrugation test were recalculated by means of weighted coefficients in accordance with the results obtained previously by subjective methodology. Six coefficients of the obtained parameters were defined using the particularly developed software for finding the smallest possible differences between the evaluated quality by means of the image analysis and with use of subjective monitoring by an experienced operator. An example of computed quality in relation to the quality of subjective assessment is showed in Figure 30.

Figure 30.

Relation between computed quality and the specification of quality by an experienced operator.

Described system has been successfully tested for two years on a production line with the accuracy of 0.1–0.3 (in the used scale 1.5–3.5) in the latest version. Final results of over a thousand measurements from a real production process show a very good potential for the on-line application [3, 4, 8, 20]. The on-line corrugation test was implemented and showed good stability of the evaluation methodology. Automatic measurement every few second is one of the many advantages. The system has been tested in a trial run.

4. Evaluation of 2D pictures with structure and defects

The explicit, objective, and automatic description of image complexity can be achieved by different methods, both statistical and using the FD. Only a few possibility examples are presented below. The analysis is done by evaluation of 2D images of surface defects—structures of the hole cracks in costume jewellery.

A digital image is represented as a matrix of pixels (or matrices for colour image, Figure 31D). Pixels can achieve different numbers, which depend on the format used for the digital images. The pixels are represented by numbers between 0 (black) and 255 (white) for the grey 8-bit palette bitmap, and the bitmap has only one matrix (brightness scale). Figure 31C shows two typical, poor quality surfaces of costume jewellery fissures. The cutting C-1 represents deep cracks, and C-2 represents a thin structure.

Figure 31.

Preparation of samples, taking photographs, software preparation.

Matlab and HarFa software [7] were used for these experimental evaluations. A methodology for analysis of the pictures was developed based on: histogram evaluation, percentage of black pixels, percentage of large defects, the FD.

In practice, the methodology of surface structure description can be divided to five steps:

• Sample preparation—the costume jewellery is cut, as the structure must be visible, Figure 31A.

• Image acquiring—cracks in costume jewellery hole. An electronic microscope was used, Figure 31B.

• Software preparation, Figure 31C (separation of image parts needed for analysis)

• Image analyses.

• Evaluation of results.

4.1. Histogram

An evolution of a bitmap structure is possible by the histogram statistical description, Figure 32. Modus, median, average, range, standard deviation, and other statistic tools can be used easily. However, they are not applicable for surface structure of costume jewellery holes.

Figure 32.

Histogram cut off on 5% level.

A potentially applicable method is the histogram cut off on 5% level. The image is described by a single number. It computes a 90% width of all the histogram pixel values from an average image value. However, this method is highly sensitive to shades that can possibly occur in the hole cracks. This easy method allows describing of all the defects, shades, cracks and structure collectively.

4.2. Thresholding

Following analyses are dealing with thresholding procedure. Thresholding is a technique that transforms grey or colour image into a binary one (black and white). For example, the binary image can be determined from the grey 8-bit palette bitmap, where black are all pixels which fulfil specific criteria, for example, 0 ≤ black ≤ 35 and all the other pixels become white (35 < white ≤ 255), Figure 33. It means, that all pixels lesser than or equal to the threshold 35 are black and greater than 35 are white. We used the threshold 35. (More than one threshold can be used or the technique for matrixes of colour images can be used too.)

Figure 33.

Binary images obtained by thresholding of grey images (8-bit palette bitmap).

If the thresholding procedure is used for all thresholds of the grey 8-bit palette image, 256 binary images are obtained. Some of the image analyses can be done for all of obtained binary images. If a single number classifies one binary image, a spectrum of dependence between single number and threshold is given (e.g. Figure 33).

4.3. Percentage of deep cracks

The method uses binary image obtained by thresholding and computes percentage of pixels with neighbouring pixels of the same value. The analysis searches black pixels (value 0) in a binary image, which have five or more neighbouring black pixels. Figure 34 shows part of boundary crack. Black pixel in Figure 34A has five neighbours, and in Figure 34B has 8. The method is suitable for a detection of relatively large and single cracks and defects that are represented by black pixels with five or more neighbouring pixels.

Figure 34.

Pixels on boundary crack.

A spectrum of dependence between percentage of black pixels with five or more neighbouring black pixels of binary images and thresholds is showed in Figure 35A. Thresholds from 10 to 50 are the most suitable for detection of large hole cracks in costume jewellery, Figure 35B. For the threshold 50, the cutting C-1 has more single cracks and defects, numerically: T_{50 C-1}=3.17% than the cutting C-2, numerically: T_{50 C-2}=0.8%.

Figure 35.

Percentage of black pixels with five or more neighbouring black pixels spectrum.

4.4. Box dimension

The method again uses binary image obtained by thresholding. The box counting method [5, 7, 9, 10] estimates FD of a structure in an image. This estimation is a single number called the box dimension D_B. The box counting method works by laying meshes of different sizes r and then counting numbers of boxes N needed to cover a binary image (Figure 36A) completely (Figure 36B, C). The power law allows us to determine a number of boxes N(r) necessary to cover the structure:

Figure 36.

Box counting method.

$N\left(r\right)=const.\cdot r^{-D_B}$E12

where D_B is the box dimension. A relation between log₂r and log₂N(r) is known as Richardson-Mandelbrot plot (Figure 36D). The box dimension can be determined by slope s of the regression line in Figure 36D:

$s=D_B=-\frac{\text{Δ}\log N(r)}{\text{Δ}\log r}$E13

The software HarFa [7] is used for the analysis and software tools developed in Matlab make data evaluation. The box dimension is multiplied by 1000 for a better confrontation.

The fractal spectrum of the cuttings C-1 and C-2 is shown in Figure 37. The box dimensions over threshold 150 are similar, as a shadow influence is significant over this value. Results of analysis for threshold 120 are:D_{B_C-1}=1429,6 (C-1) a D_{B_C-2}=1562,4 (C-2), where the higher value represents more complex structure in the image. The cutting C-2 is more structured than the C-1 and box dimension quantifies the structures.

Figure 37.

Fractal spectrum.

4.5. Example: detection of surface structure

The structures of the hole cracks in costume jewellery were analysed. Methodology of the image analysis, mentioned in the previous chapter, analysis was developed based on: histogram evaluation, percentage of black pixels, percentage of large defects, the FD. Only the last two analyses are suitable for describing these kinds of structures.

The percentage of large defects is suitable for the detection of individual, relatively large cracks and defects, Figure 38.

Figure 38.

Computation of the percentage of pixels with neighbouring pixels of the same value—percentage of large defects.

For this reason, many analyses for the FD estimation can be used. However, the most suitable for these kinds of data and structures is the box counting method (in Figure 39). For research purposes, the dimension was multiplied by 1000.Although an automatic classification of glass defects is from the laboratory, presented example shows the possibilities of this method, its sensitivity to the different shapes of glass defects.

Figure 39.

Description of glass defects using the box dimension D_B 1000 (dimension multiplied by 1000), percentage of large defects and automatic classification.

4.6. Example: defects classification

Automatic defect classification is another example of the FD application. A classification mark of obtained images can be done automatically and can be used for monitoring and control of production system. Software HarFa [7] was used for these experimental evaluations. The developed methodology was also based on: histogram evaluation, percentage of black pixels, percentage of large defects, the estimating of FD. On the base of our research, the last two analyses were chosen as suitable for describing these kinds of structures, Figure 40.

Figure 40.

Description of glass defects using the box dimension, percentage of large defects and automatic classification.

An automatic classification of glass defects can most easily be realized by using the defined boundaries for the analysis results. The simulated classification is on the right side of Figure 40. Presented example shows possibilities of this method and its sensitivity to the different shapes of glass defects. However, it does not solve problems with on-line obtaining of the images.

5. Conclusion

The examples of used methods show possibilities of application in industry and production laboratories. Structured surface, complex time series (topologically one dimensional signals), difficulty to describe dividing curves are much more common than it can be expected. The FD is then important for an objective quantification of complexity and should be use as an additional tool for the overall analysis.

The FD is widely used in science; nevertheless, industrial applications are rather rare. Data analysis using the FD has great potential in combination with statistical and other measurements in industry. However, FD cannot substitute standard analysis tools.

Although the text is focusing to a glass industry, emphasis was put on general application possibilities, where obtained knowledge, methodology and principles for product and industrial data evaluation more widely in industry might be used [11].

Acknowledgments

The results of this project LO1201 were obtained with co-funding from the Ministry of Education, Youth and Sports as part of the targeted support from the programme entitled “Národní program udržitelnosti I”.