TOWARDS A DIGITAL TWIN OF A STORAGE TANK USING LASER SCAN DATA

This paper proposes a methodology to automatically extract components of an oil storage tank from terrestrial laser scanning (TLS) point clouds, and subsequently to create a three-dimensional (3D) solid model of the tank for numerical simulation. The proposed method is integrated into a smart analysis layer of a digital twin platform consisting of three main layers: (1) smart analysis, (2) data storage, and (3) visualisation and user interaction. In this proposed method, primary components of the tank were automatically extracted in a consecutive order from a shell wall to roof and floor. Voxel-based RANSAC is employed to extract voxels containing point clouds of the shell wall, while a valley-peak-valley pattern based on kernel density estimation is implemented to remove outlier points within voxels representing to the shell wall and re-extract data points within voxels adjoined to the shell wall. Moreover, octreebased region growing is employed to extract a roof and floor from remaining point clouds. An experimental showed that the proposed framework successfully extracted all primary components of the tank and created a 3D solid model of the tank automatically. Resulting point clouds of the shell wall were directly used for estimating deformation and a 3D solid model was imported into finite element analysis (FEA) software to assess the tank in terms of stress-strain. The demonstration shows that TLS point clouds can play an important role in developing the digital twin of the oil storage tank.


INTRODUCTION
Laser scanning, also known as Light Detection and Ranging (LiDAR), has been used to capture three-dimensional (3D) topographic data of visible surfaces of objects quickly and accurately. Particularly, a terrestrial laser scanner (TLS) captures great details of the surfaces with millimetre accuracy, which are to be widely used for creating three-dimensional (3D) geometric model, and structural engineering, for example structural inspection , deformation modelling (Truong-Hong and Lindenbergh, 2019), and finite element mesh (Kassotakis et al., 2020).
In measuring deformation of a tank shell, Pukanská et al. (2014) manually extract data points of the tank shell and used point clouds of a top and bottom sections with the section thickness of 0.1m to fit a cylinder representing to the shell. Huadong et al. (2017) sliced point cloud of the shell wall based on known elevations of courses, and a least square method is used to fit the cylinder for the points of each layer. In attempt to automatically extract the points of the shell, Truong-Hong et al. (2020) used voxel-based RANSAC to extract the shell wall for a TLS point cloud capturing from exterior. Moreover, to assess the tank in terms of structural analysis, Wu et al. (2020) used an as-design model in finite element analysis (FEA) to simulate the temperature distribution due to fire. However, the use of the asdesign model in FEA can affect to results because there is discrepancy between the as-design and as-is models (Truong-Hong and Laefer, 2013).
As a large amount of crude oil stores in oil storage tanks, which can be up to ten thousand of cubic, any accident occurred can trigger catastrophic consequences for both property and environment. As such, inspection and monitoring radically change of the tank would be kept an important role to prevent any risk. A digital twin would be a great platform for engineers able to keep track of the tank real time quantities through deformations, stress and strain. However, converting a physical model to a digital one of the tank is still manual work with high cost and time consuming. Thus, this paper presents an overall framework of a digital twin for oil storage tanks. This study mainly focuses on extracting primary components (a shell wall, roof, and floor) of the tank from a TLS point cloud automatically. Subsequently, point clouds of these components are used to estimate deformation though their point clouds, and automatically create a 3D solid model directly importing into FEA to assess the tank via stress-strain.

PROPOSED METHOD
A goal of the study is to develop a digital twin framework of oil storage tanks, which consists of two main parts: Part 1tank component extraction and Part 2 -a digital twin (Figure 1). Part 1 focusses on extracting TLS point clouds of components of the, in which primary components (a shell wall, roof, and floor) are objective, while other components (e.g., columns, girders, and rafters) are parts of future work. Part II proposes a digital twin framework consisting of three layers: (1) Smart analysis, (2) Data Storage, and (3) Visualisation and User Interaction. In the digital twin platform, this paper focusses on the smart analysis layer to explore ability to use point clouds for 3D geometric model reconstruction, deformation analysis and FEA in an automated manner.
For Part I, a point cloud P = (pi = (pi.x, pi.y, pi.z 1,Nv]) by using voxel grid with a predefined voxel size (vs0) (Step 1.1) ( Figure 2b). Subsequently, a voxel is classified as "full" (Vf), if the voxel occupies the number of points larger than a predefined minimum number of points (vmin_pts = 10); otherwise, it classified as "empty" (Ve). In Step 1.2, points pi  vi Vf (full voxels) are assumed to represent a plane, as such principal component analysis (PCA) (Hoppe et al., 1992) is employed to estimate a local plane describes as vi.s = (vi.p0, vi.n), where vi.p0 is a centroid of the points pi, and vi.n is the normal vector of a fitting plan. Next, the voxels Vf are respectively classified as vertical (Vv) and non-vertical voxels Vh (Equation 1), which is based on deviation between the normal vector of the voxel and an unit vector of the axis oz (Figure 2c).
where nz = a unit vector of oz axis th = the predefined angle criterion to determine the voxels contain points of vertical surfaces  Step 2 are extracted data points of the tank shell in two sub-steps. In Step 2.1, RANSAC based voxels proposed by Truong-Hong et al. (2020) was employed to estimate the cylinder representing the shell wall from Vv, in which the cylinder (Cv) parameters consist of a center cc = (xc, yc, zc), a radius rc, a directional vector of an axis tc = (tx, ty, tz). Next, two consecutive filtering steps called normal and radius filters are used to extract the voxels occupied the points of the shell (Equation 2 and 3) (Figure 2d and 2e). The first filter is based the principle is that the point/voxel on the cylinder have a normal perpendicular to the cylinder surface, while the second filter is based on that the voxels on the cylinder surface have nearly the same radius.
The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume XLVI-4/W4-2021 16th 3D GeoInfo Conference 2021, 11-14 October 2021, New York City, USA where tsvi = the vector connected the voxel center (vi.p0) and its projection on the cylinder axis defined by a center cc and tc n = an angle criterion for a deviation between the normal of voxels and a normal vector of the cylinder vni.r = the distance from vni .p0 to the cylinder axis v.r1 and v.r2 = a lower and upper bound of the largest peak of KDE generated from vni.r, which is determined by using valley-peak-valley pattern (Truong-Hong and Lindenbergh, 2021).
As vsi may contain points of components (e.g., a roof or pipes) attached to the shell wall, and points of the shell wall may locate with other voxels (vnsj  V\Vs) that are not recognized as vsi, Step 2.2 filters those points in two sub-steps. First, for each vsi, a radial coordinate pir, the distance from the point to the cylinder axis, is computed, and the point pi  vsi belongs to the shell wall points (Ps) (Figure 2f p.r1 and p.r2 = a lower and upper bound of the largest peak of KDE generated from pi.r Second, to extract points of the shell wall located within other voxels vnsj, the algorithm starts from the voxel vnsj to search adjoined voxels as vsi  Vs, and pi  vnsj is considered as the point of the shell wall points (Figure 2g A octree-based region growing (Vo et al., 2015) is employed to segment the voxel vi, in which pi  Ps within the voxel vi is no longer consideration. To improve efficiency, only points Prp close to the top and bottom of the shell wall are respectively used as input data set for extracting the roof (Pr) and floor (Pf), which is the largest segment (Step 3.1 and 3.2) (Figure 2j).
In Part 2, to create a 3D model of the shell wall, a point cloud Ps is linearly sliced along a cylinder direction. At each cross-section, the points within the buffer (sb) are retrieved and a least square method is employed to fit a circle through these points ( Figure  3a). Next, the sweeping technique is employed to generate the 3D model based on a set of the fitting circles and known thickness of the courses (Figure 3b). This methodology is used because: (1) to avoid high cost to fit the cylinder of the shell wall, and (2) to give a better shape of the shell wall when the shell wall is subjected to deformation. Notably, as thickness of the shell wall cannot automatically identify from laser scanning point cloud in this study, those values obtained from exiting documents (e.g., design documents) are used as input data. Similarly, this procedure is also applied to the points Pr of the roof and Pf of the floor to create their 3D models, in which the cross-section is created along the radial direction.
Next, to analyse deformation of the shell wall, best fit vertical deformation and roundness analysis according to requirements of API (American Petroleum Institute, 2014) are computed using Equation 5.
bestfit/roundness = d(pi, Cv) = pi.r -rc (5) where Cv(cc, rc, tc) = the cylinder parameters consisting of a center cc, a radius rc, and an axis direction tc pi.r = a distance from the point pi to the cylinder axis For the best fit vertical deformation analysis, as the shell wall is assumed as the perfect vertical cylinder, cc and rc are derived from the fitting circle based on a projection of Ps onto the xy plane, and tc = [0, 0, 1]. For roundness, the similar procedure in Step 2.1 is used to determine the cylinder parameters. Figure 4 shows deformation from best fit vertical and roundness analyses.  Additionally, the 3D geometric model is also directly imported to FEA software to analyse behaviour of the tanks under different load scenarios to determine behaviour.

EXPERIMENTAL TEST, RESULTS AND DISCUSSIONS
To demonstrate the proposed framework, a storage tank with a diameter of 16.4592 m and height of 11.44m, located at British Columbia, Canada is selected as a case study. The tank was scanned from an interior by using Leica Scanstation P40 (Leica Geosystems, 2021). After registration, the point clouds were down-sampled to a sampling step of 5mm and 22.9 million points were exported for next steps (Figure 5a).
To extract components (a shell wall, roof and floor) of a storage tank, following input parameters were: (1) for generating a voxel grid: a voxel size vs0 = 0.2m and a predefined minimum number of points vmin_pts = 10; (2) for classifying voxels as vertical or horizontal: an angle threshold th = 45 degrees; (3) for estimating cylinder parameters using RANSAC-based voxel (Truong-Hong et al., 2020): the number of voxels to estimate the model vmodel = 500, a tolerance to identify inlier voxels max = 0.01m, and the angle criterion (n = 5 degrees) for a deviation between the normal vector of a voxel and normal of the cylinder. Moreover, for filtering outlier points of the shell wall and obtain points of the shell located within vnsj ( Step 2), the bandwidth (bw) of 0.01m is used for generating KDE. Next, for input parameters for octreebased region growing (Vo et al., 2015) are respectively an angle, distance, and residual threshold v = 5 degrees, dv = 0.01m, and rv = 0.01m. The proposed method shows that all desired components (a shell wall, roof, and floor) of a tank are extracted successfully in a term of object (Figure 5b and 5c). To assess deformation of the shell wall, the point cloud (Ps) was used to estimate the cylinder parameters by using the procedure proposed by Truong-Hong et al. (2020), in which the points were used instead of the voxels as in this proposed method. Parameters of the cylinder are (1) for the best vertical deformation cc = [-5.878, 2.637, 0.000], rc = 8.220, tc = [0, 0, 1], and for roundness cc = [-5.878, 2.639, 5.473], rc = 8.220, tc = [0.002, 0.001, 1.000].
Resulting radius of the shell wall for both analyses are nearly the same and differ from the radius from design documents no more than 9.6 mm (8.2296 m vs. 8.220 m). The difference can cause by either errors due to data acquisition and point cloud processing or damage of the shell wall or both. Moreover, based on the cylinder parameters, deformations of the shell wall at each data point were shown in Figure 5d and 5e. The shell wall was subjected to deformation, in which the minimum and maximum deformation were a range from -66.0 mm to 52.0mm for the best fit vertical analysis, and from -66.0 mm to 47.0 mm for the roundness. Moreover, the roundness analysis shows that the shell wall was inclined 0.14 degree compared to the vertical direction (tc = [0.002, 0.001, 1.000] vs. nz = [0.0, 0.0, 1.0]).
Additionally, the point clouds of the shell wall and floor were used to create a 3D model of the tank to demonstrate ability to use a point cloud for finite element analysis. As thickness of the shell wall varies along a vertical direction but is a constant for each course, the point cloud (Ps) of the shell wall were sliced along the cylinder direction. In this study, each course was divided into 4 cross-sections with an equal interval. Notably, as thickness of the shell wall and bounds of courses of the shell wall cannot automatically identify from laser scanning point cloud in this study, those values obtained from exiting documents (e.g., design drawings) are used as input data (Table 1). At each crosssection, the buffer sb = 0.01m (two times of a sampling step) were used to extract the points of the cross-section, and then a fitting circle was created (Figure 6a and Table 2). Subsequently, the sweeping technique was employed to create a 3D solid model of the shell wall. Similarly, the cross-sections along a radial direction of the floor were created, in which the interval space of 1.0m was used, and then the solid model of the floor was created, in which the thickness of 7.938 mm was used. Resulting solid models of the shell wall and floor were shown in Figure 6b.   (Table 3). This showed that the cross-section of the shell wall is subjected to minor damage and low proportion of outlier points within a sub-dataset of the cross-section. Moreover, comparing to the radius derived from design documents, radii of fitting circles differ from ones of the design document no more than 51.6mm (an average difference of 14.0mm), while RMSE is about 18.2mm.
The solid model of the tank was subsequently imported into to ANSYS software (Ansys, 2021a). A 3D element SOLID185 having 8 nodes associated with 3 degrees of freedom, and isotropic behaviour was used to model the shell wall and floor (Ansys, 2021b). The FEM of the tank consists of 49,943 nodes and 38,047 elements. The tank components (a shell wall, floor, and roof) were made from steel with mechanical properties including Young's Modulus of 210GPa and Poisson ration v = 0.3. Additionally, the tank floor was seat on gravel with Young's Modulus of 20MPa. As the goal of this analysis is to demonstrate ability of using the solid model generating directly from point cloud for FEA, the tank is assumed that crude oil filled a full tank. As such, the tank is subjected to hydraulic pressure defined the mass density (crudeoil = 875.5 kg/m 3 ) and the level of crude from the shell wall top, and the boundary condition between the floor and foundation was modelled as elastic support. Results of total deformation, equivalent stress and strain were shown in Figure 7, in which maximum total displacements, equivalent stress and strain are respectively 5.05mm, 114 MPa and 1.6x10 -3 m/m.

CONCLUSIONS
This paper proposes a digital twin framework for oil storage tanks, in which the geometric model of the tank was created from TLS point clouds. The digital twin consists of two main parts: Part 1: extract point clouds of components of the tank, and Part 2: simulate and visualize information of the tank. This study focuses on (1) extracting a shell wall, roof, and floor, (2) analysing deformation of the shell wall, and (3) demonstrate behavior of the tank via finite element analysis, in which the solid model of the tank was automatically created from its point clouds.
To extract the point cloud of the shell wall, the proposed method called RANSAC-based voxel method developed by Truong-hong et al. (2020), originally developed to extract data points of the shell wall from the data set capturing from an exterior side of the tank, was refined. This study introduced a methodology to remove outlier points using a valley-peak-valley pattern of KDE, which allows to improve quality of extraction of shell wall. Moreover, octree-based region growing was employed to extract the roof and floor. The methodology to extract components in sequence order can reduce the complexity of the data set, and the methods were based on the voxels created from sub-sets allowing to performance of the proposed method.
The proposed framework successfully extracted point clouds of the tank components, which are subsequently used to measure deformation, and automatically create a 3D model. Best fit vertical and roundness analyses used to assess deformation of the shell wall according to API requirements were demonstrated. Roundness analysis showed that the shell wall was subjected to deformation with maximum inward and outward deformation of 66.0mm and 47.0mm respectively, and the shell wall was not perfectly vertical as a design stage, in which a current inclined angle of the shell wall is 0.14 degree. Moreover, the deformation can explain a difference of a cylinder radius created from a point cloud of the shell wall and design documents about 9.6mm. A radius of a fitting circle based on point cloud of each crosssection of the shell wall is average 14.0mm differing from the radius of the designed shell wall. Finally, finite element analysis with full filled crude oil demonstrated ability to use the solid model from the point cloud to simulate behavior of the tank.
The success of this study can offer an alternative solution to assess deformation damage of an oil tank and create a 3D solid model for simulating behaviour of the tank in automated manner. The proposed method can be integrated into a digital twin platform allowing to monitoring status of the tank nearly realtime. However, other components (e.g., columns, girders, and rafters) of the tank must be extracted and assessed. Moreover, the high accuracy of behavior of the tank can be achieved if all components are included in the FE model. In future work, the methodologies for storing different data types such as point clouds, geometric models and resulting analyses and retrieving them for visualisation and user interaction will be investigated. Moreover, impact of point cloud accuracy to results of deformation measurement and of simulation will be evaluated.