Abstract
Systems in service continue to degrade with the passage of time. Pipelines are among the most common systems that wear away with usage. For public safety, it is of utmost importance to monitor pipelines. Magnetic flux leakage (MFL) testing is a widely used nondestructive evaluation (NDE) technique for defect detections within the pipelines, particularly those composed of ferromagnetic materials. Pipeline inspection gauge (PIG) procedure based on line scans can collect accurate MFL readings for defect detection. However, in real world, applications involving large pipe sectors such as extensive scanning techniques are extremely time consuming and costly. In this article, we develop a fast and cheap methodology that does not need MFL readings at all the points used in traditional PIG procedures but conducts defect detection with similar accuracy. We consider an under-sampling based scheme that collects MFL at uniformly chosen random scan points over large lattices instead of extensive PIG scans over all lattice points. On the basis of readings from the chosen random scan points, we use kriging to reconstruct MFL readings. Thereafter, we use thresholding-based segmentation on the reconstructed data for detecting defective areas. We demonstrate the applicability of our methodology on synthetic data generated using finite element models and on MFL data collected via laboratory experiments. In these experiments, spanning a wide range of defect types, our proposed novel MFL-based NDE methodology is witnessed to have operating characteristics within the acceptable threshold of PIG-based traditional methods and thus provide an extremely cost-effective, fast procedure with competing error rates.
1 Introduction
Monitoring the condition and performing effective diagnosis of the defects within the pipeline is a necessity as the structural integrity of the pipeline decreases with time. Many oil and gas pipeline failures have led to fatalities and significant loss of properties in the recent past; hence developing inspection techniques to access their conditions is of considerable importance [1–3]. Since the 1960s, magnetic flux leakage (MFL) serves as one of the most widely used nondestructive evaluation (NDE) techniques for in-service inspection of pipelines as the pipeline materials are mostly ferromagnetic [4–6]. The working principle of MFL is as follows: ferromagnetic pipe wall is magnetized close to saturation by the aid of a permanent magnet or a coil wound on a ferromagnetic yoke; the presence of defects decreases the wall thickness; and so, magnetic flux density and reluctance are increased in the vicinity of defects [7]. A higher fraction of the magnetic flux will thus leak from the pipe walls near defects into the air. Advanced MFL signal-processing techniques that can detect defects with high accuracy across different scenarios have been developed in the recent literature [8,9]. While these MFL-based methods have been very successful to detect defects in noisy situations [9] as well as for dynamic tracking of defect sizes using transfer learning on sequential scans [10], they cannot be used in monitoring large pipelines as they involve near-continuous line scans, which is extremely time consuming. For monitoring massive pipelines, here we develop a fast-approximate algorithm that uses a minuscule fraction of data but provides similar operational performance as the aforementioned data-intensive MFL methods.
In this article, an under-sampling-based data processing scheme via kriging [11] is proposed. The procedure works without scanning the entire massive pipelines and only using MFL readings from a very few random scan points in these pipelines, thereby saving humungous time and cost. The key contributions of this article are as follows:
By using kriging, we conduct spatial interpolation by using the limited set of sampled data points to estimate the magnetic flux over the concerned continuous pipe sector.
Scan points in the defective areas have very different MFL features than those from nondefect regions. Here, the method considers the spatial autocorrelation among the sampled data points, and in this process, it allows tracking and characterization of defects from much fewer readings than those permitted in pipeline inspection gauge (PIG) based on continuous scanning.
The present work involves developing a three-dimensional (3D) finite element (FEM) model in comsol to investigate the performance of the proposed technique in detecting the surface and subsurface defects of various shapes and sizes. See Sec. 2.1 for details.
The efficacy of our proposed kriging-based fast and approximate method is tested against a wide range of defect categories such as single, multiple, and interacting defects. Thresholding is used to segment defective and nondefective areas in pipelines based on MFL readings from (a) exhaustive line scans and (b) kriging interpolated predictions based on significant smaller sample of scan points. Encouraging results (see Table 2 in Sec. 4) are obtained, which shows the predicted defective areas by our proposed method to be in close resemblance with the ground truth.
The results produced on the simulated data were also validated experimentally (see Table 3 in Sec. 4) by developing a setup in the laboratory (see Sec. 2.2) that enabled collecting MFL readings akin to real defects.
Term | Length | Width | Height | Material |
---|---|---|---|---|
Yoke | 400 | 50 | 40 | Ferro-nickel alloy |
Magnet | 30 | 80 | 40 | NdFeB |
Brush | 30 | 50 | 40 | Ferro-nickel alloy |
Specimen | 600 | 400 | 10 | X52 like iron |
Term | Length | Width | Height | Material |
---|---|---|---|---|
Yoke | 400 | 50 | 40 | Ferro-nickel alloy |
Magnet | 30 | 80 | 40 | NdFeB |
Brush | 30 | 50 | 40 | Ferro-nickel alloy |
Specimen | 600 | 400 | 10 | X52 like iron |
Defect type | Performance measure | Subsampling rate (%) | |||
---|---|---|---|---|---|
5/1 | 10/1 | 20/1 | 30/1 | ||
Case I | R2 | 92.41 | 65.67 | 27.07 | 13.75 |
Coverage | 83.93 | 51.79 | 25.00 | 7.14 | |
Exceedance | 0.00 | 0.25 | 0.19 | 0.06 | |
Case II | R2 | 88.67 | 65.48 | 47.40 | 9.84 |
Coverage | 73.02 | 69.84 | 19.05 | 1.59 | |
Exceedance | 0.38 | 0.50 | 0.13 | 0.00 | |
Case III | R2 | 93.55 | 76.41 | 36.82 | 9.85 |
Coverage | 96.67 | 82.33 | 73.67 | 6.67 | |
Exceedance | 0.38 | 1.94 | 5.38 | 0.00 | |
Case IV | R2 | 90.17 | 77.80 | 52.56 | 20.00 |
Coverage | 90.91 | 71.16 | 68.97 | 6.90 | |
Exceedance | 0.50 | 1.38 | 5.38 | 0.44 | |
Case V | R2 | 96.11 | 88.14 | 70.16 | 27.56 |
Coverage | 89.83 | 84.88 | 71.51 | 8.43 | |
Exceedance | 0.56 | 1.00 | 3.94 | 0.50 |
Defect type | Performance measure | Subsampling rate (%) | |||
---|---|---|---|---|---|
5/1 | 10/1 | 20/1 | 30/1 | ||
Case I | R2 | 92.41 | 65.67 | 27.07 | 13.75 |
Coverage | 83.93 | 51.79 | 25.00 | 7.14 | |
Exceedance | 0.00 | 0.25 | 0.19 | 0.06 | |
Case II | R2 | 88.67 | 65.48 | 47.40 | 9.84 |
Coverage | 73.02 | 69.84 | 19.05 | 1.59 | |
Exceedance | 0.38 | 0.50 | 0.13 | 0.00 | |
Case III | R2 | 93.55 | 76.41 | 36.82 | 9.85 |
Coverage | 96.67 | 82.33 | 73.67 | 6.67 | |
Exceedance | 0.38 | 1.94 | 5.38 | 0.00 | |
Case IV | R2 | 90.17 | 77.80 | 52.56 | 20.00 |
Coverage | 90.91 | 71.16 | 68.97 | 6.90 | |
Exceedance | 0.50 | 1.38 | 5.38 | 0.44 | |
Case V | R2 | 96.11 | 88.14 | 70.16 | 27.56 |
Coverage | 89.83 | 84.88 | 71.51 | 8.43 | |
Exceedance | 0.56 | 1.00 | 3.94 | 0.50 |
Defect type | Performance measure | Subsampling rate (%) | |||
---|---|---|---|---|---|
5/2 | 5/1 | 15/2 | 10/1 | ||
Rectangle | R2 | 99.42 | 94.6 | 85.01 | 70.95 |
Coverage | 99.34 | 95.8 | 81.24 | 70.42 | |
Exceedance | 0.32 | 1.52 | 0.52 | 3.68 | |
Sphere | R2 | 99.51 | 98.83 | 96.70 | 80.47 |
Coverage | 100.00 | 100.00 | 99.70 | 80.78 | |
Exceedance | 0.40 | 1.76 | 3.24 | 0.72 |
Defect type | Performance measure | Subsampling rate (%) | |||
---|---|---|---|---|---|
5/2 | 5/1 | 15/2 | 10/1 | ||
Rectangle | R2 | 99.42 | 94.6 | 85.01 | 70.95 |
Coverage | 99.34 | 95.8 | 81.24 | 70.42 | |
Exceedance | 0.32 | 1.52 | 0.52 | 3.68 | |
Sphere | R2 | 99.51 | 98.83 | 96.70 | 80.47 |
Coverage | 100.00 | 100.00 | 99.70 | 80.78 | |
Exceedance | 0.40 | 1.76 | 3.24 | 0.72 |
In Section 2, we describe simulation procedures and lab experimental setups for acquiring MFL data akin to real-life scenarios. We consider various types of interacting defects in pipe sectors, which possess more threat than the singular defects. In Sec. 3, we describe our the kriging-based fast NDE procedure and its advantages. In Sec. 4, we present the applicability of our proposed methodology. Our results in Sec. 4 show that by using only a very small fraction of MFL readings from the datasets of Sec. 2, our proposed method can attain accuracy as good as extensive scan-based traditional PIG methods. In Sec. 5, we end with a discussion on the prescribed method.
2 Magnetic Flux Data Generation
For both synthetic and experimental data in this section, we generate MFL readings for all the scan points in a metallic pipe sector. To generate synthetic data, we design a 3D finite element MFL model in the comsol multiphysics modeling software. The synthetic data generation procedure is also described. Thereafter, we present the experimental procedures.
2.1 Simulation Designs.
Our simulation data generation process is constructed as a magneto-static problem governed by the conventional Maxwell’s equation. Figure 1 shows the different components of a 3D MFL model in comsol along with a defective sample with five interacting defects of varied volumes. The different dimensions and the material properties needed to construct the model is presented in Table 1. The benchmark as discussed in Ref. [10] is used here to construct the model with the lift-off parameter being set at 2 mm.
As shown in Fig. 1(b), a rectangular surface consisting of multiple rectangular notch defects has been considered. Different scenarios like varying the size of notches, interaction among the notches, and evolution of new defects are considered. Parametric sweep in x and y direction is conducted as shown in Fig. 2, and the magnetic flux in the axial direction Bx is collected in the form of matrices. Near the defects, extremely fine triangular mesh is used, whereas on the rest of the surface, tetrahedral fine mesh is employed.
We use a rectangular grid of length 8 cm and breadth 5 cm. We consider heterogeneous spacing along the x and the y directions in the grid. The spacing between grid points is 0.5 mm along the x-axis and 0.2 mm along the y-axis. Consequently, we have 160 × 250 = 40,000 grid points. We initially consider four interacting defects of equal sizes. Figure 2(a) shows the XY plot of a two-dimensional (2D) fine rectangular grid containing the location of these four defects. The design is symmetric and balanced in defect sizes. Later we consider situations where another defect is added to a new location on the rectangular lattice and some of the existing defects evolve to grow larger in sizes. Figure 2(b) shows the schematic for this case with the upper ones increasing in size and a new defect cropping up on later inspection.
We generate MFL readings pertaining to five different designs. With the design in Fig. 2(a) being the baseline, we allow the defects in the subsequent four cases to either increase in size or remain constant. The simulations reflect the perturbation in MFL readings due to not only the increase in defect sizes but also due to interactions in flux leakages from the neighboring defects. The 2D fluctuations of the magnetic fields along the axial direction Bx due to spatial movement along x and y axes are observed in comsol from these designs. Figure 3 plots these MFL readings. From the figure, it is evident that in the vicinity of the defects, the perturbation in the magnetic field is significantly higher than those from the nondefective points.
The five different designs considered here are as follows:
Case I: All the defects are squares with the length of each side being 5 mm. Their locations in the rectangular lattice of length 2.5 cm and width 1.6 cm are given by the layout in Fig. 2(a).
Case II: The upper two defects in case I have increased to squares of length 9 mm, whereas the bottom three defects are squares of length 5 mm each. A new defect has evolved in case II compared to case I. The layout is given in Fig. 2(b).
Case III: The upper defects have increased in size to squares whose length of each side is now 14 mm.
Case IV: The upper defects have increased to squares with 18 mm sides. The MFL reading plot in Fig. 3 shows increasing interaction among the upper defects.
Case V: The upper defects have merged into a single defect on further increasing their dimension to 20 mm sides each. The lower defects remain constant in size.
Figure 3 shows the plots of the MFL data generated in comsol from the aforementioned five different designs.
2.2 Experimental Setup.
The working principle of MFL-based NDE methods is based on the fact that there will be an increased magnetic flux density in the vicinity of the defects [6], which can be measured using a hall effect sensor (static and dynamic fields) or coil (dynamic fields) [7]. Following this principle, we design our experimental setup using the following ingredients:
A permanent magnet-based MFL probe.
A scanning robot arm to move the probe along the sample.
A direct current (DC) power supply.
A data acquisition system (DAS).
Other associated units.
The image of the setup and its associated MFL probe is shown in Figs. 4(a) and 4(b), respectively. As shown in Fig. 4(b), the permanent magnet assembly includes a magnetic circuit consisting of two permanent magnets (NdFeB) and the magnetic field sensor to sense the leaked magnetic field in the presence of defects. The robotic scanning arm, which holds the MFL probe, is used to control the movement of the permanent magnet assembly. We use an analog giant magnetoresistance (GMR) magnetometer sensor (AAH002-02E) manufactured by NVE Corporation in our sensor setup. Figure 5 shows the schematic of the entire experimental setup, the pin configuration of the GMR probe, and the data acquisition system. The GMR sensor contains four resistors in the form of Wheatstone bridge configuration [12,13]. We chose AAH series of sensors as it has high sensitivity for low-field sensing and excellent temperature stability. Also, the small size of the sensor makes it very convenient for mounting in the constructed MFL probe.
The axis of sensitivity of the sensor is parallel to the surface of the test material. We record the changes in the axial component of the magnetic flux (Bx) based on the output voltage of the GMR sensor. When the MFL probe is far from the ferromagnetic sample, then the output voltage of the GMR sensor is constant. In the presence of defects, the magnetization changes that subsequently changes (Bx). This produces a resistance change in the GMR sensor, thereby altering the output voltage. This voltage is subsequently analyzed by our proposed algorithm. A constant DC voltage of 1.98 V is given as input to our experiment. This low-field range and high resolution of the GMR sensor make it ideally suited to measure the residual fields. The pivotal ingredient in our DAS was the National Instrument Data Acquisition card (PCIe-6341), which samples and digitizes the data using an imaging routine and the output is plotted on a computer. The sensor was interfaced with the data acquisition card via digital multimeter to allow simultaneous data acquisition.
We use square steel samples with 50 cm sides. We consider two different experiments one with a circular and the other with a rectangular defect. Both the defects were placed in the center of the steel samples. The rectangular defect is of length 2.5 cm and breadth 1.2 cm, whereas the spherical defect has a diameter 0.7 cm. Our designed MFL probe recorded flux readings at scan points 1 mm apart producing a 50 × 50 square grid of MFL records. Their plots are shown in Fig. 6 (rectangular defect) and 7 (spherical defect).
3 Kriging-Based Proposed Methodology
3.1 Magnetic Flux Leakage Prediction.
Kriging is a method of spatial interpolation that originated in the field of mining [11]. Precision provided by the state-of-the-art PIG techniques comes at a price of computational cost. In large-scale NDE applications, this cost becomes crucial and can be prohibitive for timely damage control. Thus, it is extremely important to consider cheap surrogate methods that reduce the computational cost, while maintaining the required precision.
We use Kriging, which is a spatial interpolation method to develop scalable NDE techniques for defect detection in large pipelines. Kriging uses a limited set of under-sampled data points to estimate the variable over a continuous spatial domain. The interpolation is based on the spatial arrangement of the empirical observations, rather than on a presumed model [14]. Thus, kriging generates estimates of uncertainty surrounding each interpolation. Kriging predictor is a combination of linear predictor and exact interpolator, thereby the value obtained by kriging for any actual sampled location will be equal to the observed value at that point, and the interpolated values will be the best linear unbiased predictor.
Magnetic flux distribution depends on the location of the coordinate points. As there exists a strong spatial correlation among the MFL from neighbor points on the lattice, kriging can provide a good prediction for the unobserved points by leveraging the correlations with their neighboring observed scan points. Kriging is basically a two-step process, where, in the first step, the spatial covariance structure of sampled points is determined by fitting a variogram, which is a visual depiction of the covariance of each pair of points in sampled data. Thereafter, the weights from this variogram are used to interpolate values in the unobserved points [15].
Let G be the grid of all points used in a traditional PIG line scan or 2D scan. Let GS be the subset of points in G that are used for kriging. Let M be matrices containing standardized magnetic flux leakages readings for points in G and MS be the set of all sampled M values at GS. Let be the predicted MFL values based on our procedure that takes MS and the location of the scan points GS as inputs. The entries and mi,j in and M are the predicted and recorded MFL values based on our proposed method and PIG line scans, respectively, at the (i, j)th point in the grid G.
3.2 Subsampling Rate and Defect Detection.
The subsampling ratio adjusted for dimension is defined as follows:
Sampling rate per dimension = ρS = (|G|/|GS|)1/d, where |G| and |GS| are respectively the number of scan points in G and GS and d represent the dimension of the grid. For a fixed threshold h, we segment points based on the standardized MFL values. The readings in Figs. 6 and 7 (experimental data) from the defective points are troughs in lattice, whereas that in comsol simulation (see Fig. 3) constitute a crest. The orientation of the crest and troughs depends on the different data collection procedures. In our simulation experiments, the variation in magnetic flux signals is calculated directly and plotted. However, from the lab experiment, the changes in the output voltage from the GMR sensor are plotted. As the defect detection algorithm used in this article only requires the readings from the defective region to be different from that of the background, we are not concerned with whether the readings of defective regions were uplifted or downregulated. Using threshold h and the complete PIG scan data, we segment the points in G into defective and nondefective or trouble-free sets D and T, respectively. For most cases, h set as 0.5 works well with standardized MFL values. Keeping the same value of h and the same criterion, we segment our Kriging predictions into defective and trouble-free sets and , respectively. Ideally, we would like and to be very similar to sets D and T. We analyze their relations in the following section.
4 Results
The larger the subsampling rate ρS, the fewer is the number of points used in our proposed method. The greater ρS is from 1, the faster the MFL data collection step in our procedure will be. However, if |GS| is too small, then we will do a shoddy job in interpolating the MFL values at the unsampled points , and any subsequent inference on defect location recognition will be highly erroneous. So, we consider popular metrics that capture the operating performance of defect detection algorithms. We report the mean square error (MSE) of the MFL prediction based on kriging. We also report the coverage of defect points based on our proposed algorithm as well as its false-positive percentages. Next, we define these measures.
Lower MSE signifies that the predicted MFL values by the kriging model is closer to the exhaustive scan MFL values, in which case the R2 values will be closer to 1, signifying accurate reconstruction of MFL values by kriging at the unsampled points.
High R2, high coverage, and low exceedance values are desired. Next, we report these metrics as subsampling rate ρS is varied. We consider uniform subsampling designs throughout the article.
4.1 Results on Synthetic Data.
Consider the synthetic datasets described in the five designs of Sec. 2.1. In Table 2, we report the coverage, exceedance, and R2 performance measures of our proposed methodology as the subsampling rate ρS varies from 5 to 30. Note that in these comsol simulated datasets, |G| = 40,000, and thus, a subsampling rate of 5:1 means 40,000/52 = 1600 random locations in G were considered in GS. For ρS equal to 10, 20, and 30, we, respectively, have 400, 100, and 44 MFL samples.
From the table, we see that when ρS = 5, we get considerably close defect identification with very low false positives across all the five cases even if we consider only of the MFL readings. The R2 values of the predicted MFL values are quite good signifying that our proposed method can be used with high confidence at ρS = 5. Figure 8 shows the MFL predictions by our method for case II of Table 2. It shows the gradual deterioration in the kriging-based MFL reconstruction as the sampling rate ρS increases. The coverage rate decreases considerably when ρS = 20, and as expected, there is a complete breakdown of the method as with ρS = 30 for there are only 40 scan points.
4.2 Results on Experimental Data.
Next, we apply our proposed method for the experimental data described in Sec. 2.2. The different performance metrics are reported in Table 3, and Figs. 9 and 10 show that predicted readings for the cases with rectangular and spherical defects, respectively, as subsampling rates are increases from 2.5 to 10.
From Figs. 9 and 10 and Table 3, we see that in both the cases, our method produces coverages sufficiently close to exhaustive scans till the subsampling rates reach 5. Thereafter, the coverage decreases below the tolerance limit for the rectangular defect. In all these cases, the exceedance is very low. In particular, from the 3D MFL plots in Figs. 9 and 10, we observe that the troughs corresponding to the rectangular and spherical defects become greatly thin and lose much of their shape as subsampling rates are increased from 5 to 7.5. In Fig. 11, the standard deviation of the estimates associated at each scan point is plotted. In Table 4, we report the mean of standard deviation (MSD) across all the scan points. Note, that as we have used standardized reading values, the range of the readings is 1. To understand the relative impact of these standard deviations, we report the relative standard deviation (RSD) by dividing MSD by the average difference between the readings from defective and nondefective scan points. Higher RSD means greater confidence in correctly predicting defective scan points. In Fig. 12, the prediction surface for the two defects is also reported for different subsampling ratios. From Table 4, we observe that for subsampling rates up to 5, the MSD and RSD are well controlled. Figure 12 shows that the prediction surface is thin enough to provide accurate differentiation between the predicted readings from defective and nondefective scan points.
Defect type | Measure | Subsampling rate | |||
---|---|---|---|---|---|
5/2 | 5/1 | 15/2 | 10/1 | ||
Rectangle | MSD | 0.0166 | 0.0622 | 0.0641 | 0.1367 |
RSD | 18.5940 | 4.9744 | 4.8274 | 2.2614 | |
Sphere | MSD | 0.0041 | 0.0185 | 0.0307 | 0.0443 |
RSD | 26.7102 | 5.9282 | 3.5651 | 2.4689 |
Defect type | Measure | Subsampling rate | |||
---|---|---|---|---|---|
5/2 | 5/1 | 15/2 | 10/1 | ||
Rectangle | MSD | 0.0166 | 0.0622 | 0.0641 | 0.1367 |
RSD | 18.5940 | 4.9744 | 4.8274 | 2.2614 | |
Sphere | MSD | 0.0041 | 0.0185 | 0.0307 | 0.0443 |
RSD | 26.7102 | 5.9282 | 3.5651 | 2.4689 |
Note: MSD and RSD for different subsampling rates are reported.
5 Discussion
On the basis of the results in Sec. 4, we observe that with much fewer number of MFL readings, by using our kriging-based method, we can identify defective areas as accurately as exhaustive line scans. As such, in Sec. 4, across a wide spectrum of synthetic and experimental datasets, we found that our kriging-based method can provide more than defect identification coverages and lower than false-positive rates when it uses only of the MFL readings used in exhaustive traditional PIG approaches. Thus, our proposed method is fast, cost-effective, and highly scalable for inspecting defects in large pipelines. In future, it will be useful to study the improvements due to nonuniform sampling designs and introspect the applicability of our prescribed method to other complex NDE tasks in pipelines such as dynamic tracking of defects.
Here, we show that our proposed kriging-based detection methodology worked well in detecting millimeter-sized defects, which is popularly done by MFL-based NDE inspections [20–23]. Further work is needed to see if aberrations in magnetic flux due to smaller defects can be recovered from under-sampled signals. In this article, inspection cost is reduced by under-sampling and thereafter reconstructing the MFL signals at all scan points by using Kriging. Compared to MFL readings from extensive scans, there is always some information loss in the Kriging-based reconstructed signal. The success of our procedure rests on the fact that this information loss does not hamper the detection of the presence and location of the defects in the metallic surface. It is to be noted that the information loss will be exacerbated if the MFL values are very noisy; this can happen if there is lift-off or gap between the probe and the inspection surface. In such situations, it will be difficult to recover defect locations from under-sampled MFL readings as the signal-to-noise ratio can be very low.
Here, in lab experiment, we considered 50 cm × 50 cm metallic plates with defects. Continuous scans over large pipe surfaces will have a severe imbalance between defective and nondefective scan points. In these cases, we can apply our proposed method frame-by-frame by first identifying rectangular frames containing defects and thereafter locating the defects in frames. Most frames will not contain any defect. As our proposed method is witnessed to provide good coverage and low exceedance in defective scan points detection, we expect the low false discovery rate and high power in such frame-by-frame detection analysis.
Acknowledgment
This work is partially supported by the US Department of Transportation Grant: Improvements to Pipeline Assessment Methods and Models to Reduce Variance (Award No. 693JK1810001).
Conflict of Interest
There are no conflicts of interest.
Data Availability Statement
The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request. The authors attest that all data for this study are included in the paper. Data provided by a third party listed in Acknowledgments.