Bubble mapping: three-dimensional visualisation of gas–liquid flow regimes using electrical tomography

It is well known that electrical tomography is capable of ‘seeing’ through opaque pipe walls and flow media. However, due to the relatively low spatial resolution, electrical tomograms are ineffective on visualisation of small bubbles or sharp interfaces between large bubbles and the liquid. These limitations give rise to ambiguity in human and/or machine perception of flow dynamics from presenting bubble cluster and blurry boundaries of large bubbles in faint colour or grey. In this paper, a binary approach, called bubble mapping, is proposed to enhance the visibility of flow regime in gas–liquid pipeline flows. With the input gas concentration tomograms reconstructed by an electrical tomography, the method takes the major visual characteristics of typical flow regimes in common recognition as a prior knowledge and replaces the conventional colour mapping with a number of bubbles and their mergence with a thresholding value. With this approach, a stack of cross-sectional tomograms by electrical tomography is transformed and displayed as a collection of individual gas bubbles. The transformation is done with the help of a lookup table indexed by bubble size and an enhanced isosurface algorithm for vivid three-dimensional visualisation. This new approach has been applied to a gas-in-water flow rig and its performance is compared with photo images taken by a high-speed camera. Reasonable agreements have been obtained for common flow regimes including bubbly flow, stratified flow, plug flow, slug flow and annular flow in horizontal pipe and bubbly flow, slug flow and annular flow in vertical pipe. Results evidence the method can provide a photorealistic flow visualisation and therefore enhance the visibility of electrical tomograms for visual identification of flow regime, particularly for opaque media in industrial flows.

is in a steady state, i.e. the flow with certain flowrate of each phase is fully developed, flow regimes can be simply categorised into a few different forms depending on pipe orientation. Flow regime map can be derived to describe possible flow regimes with certain pipe orientation. Figure 1 demonstrates flow regime maps and typical flow regimes for vertical and horizontal gas-liquid pipeline flow.
The dynamics of multiphase flow are of vital importance since they significantly influence the optimal design parameters and the operational conditions of the industries containing multiphase flow [3]. Sudden changes of flow patterns in production pipe could have adverse effect on productivity, the pipe itself, and the instruments associated with the production system. Therefore, flow regime visualisation has attracted many attentions from scientists and engineers. The method of observation with a high-speed camera may be the most common one to reveal the spatial distribution of gas in water. This observation, however, is subject to the availability of the optical transparency of pipe wall and the continuous phase in fluids. Moreover, the void fraction of gas phase also affects the quality of the observation (e.g. 10% [4]). Other alternative non-intrusive method, e.g. gamma-ray tomography, is rarely applied due to its limits as its cost, installation and safety. Therefore, an effective and feasible method is highly demanded.
Process tomography has been widely applied to measure and visualise multiphase flow in the past few decades. Particularly, electrical tomography is a type of visualisation techniques based on the phase difference in electrical properties, e.g. conductivity or permittivity [3]. Compared with other tomographic techniques, electrical tomography, e.g. electrical resistance tomography (ERT) or electrical capacitance tomography (ECT), is capable of producing cross-sectional images with a relatively high temporal resolution (sub-millisecond) [5][6][7] but a relatively low spatial resolution ( e.g. around 5%, i.e. the ratio of the smallest size of identifiable objects to the diameter of the sensor) [8]. Electrical resistance tomography, as a typical electrical tomography, includes a sensor consist of a number of electrodes for injecting current into and measuring voltage drops from the boundary of process medium under the investigation, a signal conditioner and data acquisition system and an algorithm to reconstruct the conductivity distribution in the medium [9]. Generally, the local conductivity in tomogram has to be translated to an engineering data, e.g. signing gas concentration to an imaging pixel in the study, which can be converted with Maxwell relationship [10]. Algorithms can be generally categorised to one-step or multistep iterative algorithms based on the computation speed, or qualitative or quantitative algorithms based on the reconstructed data accuracy [9]. The third type of algorithms may include those based on a prior knowledge [11], statistical [12,13] or artificial intelligent [14] approaches, which presents their specific features and limits. In this study, we are focusing on the enhancement of flow visualisation from reconstructed tomograms which were reconstructed with the modified sensitivity back-projection algorithm (MSBP), a one-step algorithm [15]. It is expected the principle of method discussed in this paper can be applied to tomograms reconstructed with other algorithm.
Among the visualisation techniques for scalar fields, the method of colour mapping is the most widely applied, owing to its simplicity in computation and implementation [16]. The  basic idea behind the colour mapping is to convert different  scalar values to different colours according to a predefined  lookup table. From visualisation and perception viewpoints,  applying a lookup table should emphasise important charac-teristics, but reduce less important or irrelevant information to the minimum. However, the colour mapping-based images in electrical tomography are unable to reveal sufficient flow characteristics, e.g. a low concentration of small bubbles in a common colour (figure 2(a)), owing to the incapability of electrical tomography to report the conductivity change over a region as stated previously. The faint change in colour mapping due to the small bubbles present a challenge for human visibility. In addition, the images are strongly affected by the definition of the lookup table, which may result in ambiguity for human and machine perception (figure 2(b)), e.g. bubble boundary.
In addition to the aforementioned approaches, there are a few different approaches for the visualisation. In [4], wiremesh sensor (WMS) was used for measuring bubble size and distribution in air-water flow in a vertical vessel, without the involvement of reconstruction process, i.e. on raw data. In Prasser's method, a bubble is reconstructed by firstly locating the centre of a bubble, and then swarming the neighbouring pixels around the centre to form the bubble. However, the result may tend to be an overestimation of bubble size. When applying to electrical tomograms, this method is unable to identify small bubbles due to the relatively low spatial resolution of the tomograms. Another example is the application of raytracing technique by [17], which works by tracing a path from an imaginary eye through each pixel in a virtual screen, and calculating the colour of the object visible through it [16]. Since bubbles in electrical tomograms are not clearly presented, the ray is unable to interact with individual objects, i.e. bubbles, and hence raytracing is not directly applicable to electrical tomograms. [18] applied conventional texture-based volume rendering technique in computer graphics to achieve real-time 3D visualisation on ECT tomograms, which essentially displays 3D volume data in a 2D image by projection. Again, this method does not manage to reflect adequate flow dynamics involving small bubbles.
In this paper, a new approach is proposed to enhance the capability of electrical tomography to visualise gas-liquid pipeline flows. It is assumed that as fully-developed flow, bubble sizes proportional to local mean volume fractions, small bubbles are spherical while large bubbles are random in shape, and the growth and collapse of bubbles are ignored [1,2]. With the assumptions, a bubble lookup table is defined to replace the conventional colour lookup table. By remeshing the original dataset with temporal and spatial information into a structured space consisting of uniform cubes, the bubble lookup table can be applied to reconstruct small bubbles with a visible size but keeping the concentration the same as the original value in the cube. Large bubbles can also be located by merging a swarm of adjacent cubes with mean volume fractions beyond a carefully-selected threshold. Later, an enhanced isosurface method is utilised to identify the boundaries between large bubbles and the liquid phase.
The rest of the paper is arranged as follows. The methodology is introduced in section 2, and the visualisation results are presented in section 3. Then, conclusions are made in section 4, along with brief discussions of the benefits and limitations of the approach.

Methodology
In our study, a 16-electrode ERT is deployed. The number of independent measurements is 114, whereas the number of the unknown pixel resistance is 316. The image reconstruction  therefore is to derive the 316 unknowns from the 114 measurements. The non-linearity makes the precise derivation of the unknowns extremely challenging, even advanced iterative reconstruction algorithm is applied. Figure 3 presents the ERT concentration tomograms of small bubbles (figures 3(a) and (b)) and large bubbles (figures 3(c) and (d)) reconstructed by MSBP and sensitivity-based conjugate gradient (SCG) method [15]. The images clearly demonstrate the limitations of electrical tomography when visualising multiphase flow. That is, they are unable to identify small bubbles below a certain size or to indicate the sharp interface between the dispersed phase (e.g. gas) and the continuous phase (e.g. water). Nevertheless, despite its low spatial resolution, it has been reported, in terms of ERT, that it is capable of handling gas concentration close to 100% for a gas-water flow [19].
Conventionally, the colour of each pixel in electrical tomograms is determined by the corresponding concentration value of the pixel. For the convenience in discussion, the concentration value at the imaging pixel is dimensionlessly defined as 0% or 100%, corresponding to a pixel fully occupied by liquid or gas. Otherwise, the pixel concentration is within the range (0%, 100%), meaning the pixel is partially filled by gas. Following the same concept, a new lookup table can be established, in which a concentration value is transferred to a spherical bubble within a cell, but the concentration value is maintained by the ratio of the areas of the spherical bubble against the cell. In order to reconstruct bubbles correctly, both spatial and temporal dimensional information have to be taken into account, by means of properly defining the dimension of interrogation cell (IC). Based on the dimensional definition of IC, the original concentration space is split into a coarser regular grid, to which a new lookup table is applied to reconstruct small bubbles, while indicating the rough locations of large bubbles. Afterwards, an enhanced isosurface algorithm is utilised to identify the boundary between large bubbles and liquid. Finally, when both small bubbles and large bubbles are generated, they are displayed by the common approach in computer graphics. Figure 4 shows the steps of the new approach, and each step will be addressed in the following sections.

Data pre-processing
Due to unavoidable noise during measurement, the converted concentration distribution from conductivity distribution may contain abnormal data, e.g. negative values. Therefore, it has to be filtered firstly. The filtered distribution is then split into IC by incorporating both spatial and temporal dimensional information into the definition of the cells. In particular, spatial information needs to consider the size of small bubbles in real situation, as well as the resolution capability of the systems, while temporal information needs to incorporate actual velocity of gas phase and data acquisition speed of the employed systems into IC definition.
Given a set of real numbers R, a concentration set is defined as: The set C may need to be filtered/normalised to fit into meaningful range, i.e. (0%, 100%). When the software [20] that was employed for the reconstruction process to generate stacked concentration tomograms provides a non-shifted mean concentration, the pixel concentrations are always within the range (0%, 100%) with some internal restrictions  of the algorithm in the software, in which case filtering/normalisation is not required. Otherwise, a simple filter is defined by (2): When the filtering in the approach is required, a shift of averaged concentration with and without it may be introduced. However, the difference is actually too little to affect the quality of visualisation, and therefore the error is treated as ignorable.

Data re-meshing
In ERT, there are two important criteria, namely spatial and temporal resolution, for evaluating the quality of ERT tomograms. They refer to the precision of ERT measurement with respect to space and time. The spatial resolution of ERT is down to 5% of vessel size [8], whereas the temporal resolution varies depending on the hardware. Accordingly, for a pipe with diameter D p , the minimum size D m is D m = 0.05 × D p . On the other hand, when a gas-liquid flow is fully developed, the largest size of small bubbles can be observed and estimated. For example, when pipe diameter is within 50 mm −200 mm, the largest size is approximately less than 20 mm [21][22][23][24]. It is worth noting that the primary purpose of remeshing is to generate small bubbles missing in ERT tomograms. Together the spatial resolution of ERT and the assumption of the small bubble size, the spatial dimension of IC can therefore be defined as: s,x and Dim IC s,y are the spatial dimension Dim IC s of IC at x and y axis respectively, D s is the maximal diameter of small bubbles, and D p is the diameter of the pipe containing the flow, m and n is the spatial grid definition of the original concentration set.
It is particularly noted that the third dimension of the input tomographic dataset generally only reflects temporal information (imaging speed), rather than spatial information. It therefore needs to be transferred to spatial information by incorporating local velocity of gas phase with sufficient data acquisition speed of the tomography. For instance, suppose the local velocity of the gas phase is 1 m s −1 , and hardware collection speed is 1000 fps, the collected 1000 frames need to reflect the one-meter movement of the gas phase if the measurement is taken for one second. That is, within the time for collecting two consecutive frames, i.e. 1/1000 = 1 ms, the flow actually moves forward 1 m s −1 × 1 ms = 1 mm. Therefore, the temporal dimension of IC can be defined as: where v is the local velocity of the gas phase, and ∆t is the temporal resolution of the original concentration set.
Ideally, the velocity should contain the one for every bubble since each bubble may have different velocity. For example, concentration, velocity field, and other information are usually utilised together to for the derivation and evo lution of bubble size and distribution in multiphase flow simulation [25]. Unfortunately, velocity field is hardly available in tomography-based flow characterisation. Instead, statisticsbased methods, e.g. cross correlation [26] are applied to derive the axial velocity as an approximation of the individual velocity of each bubble, with the assumption that the velocities at other directions are ignored. Given an ERT system with two parallel-installed electrode arrays along pipeline, gas velocity can be derived with a direct cross correlation method using the following equation: (5) where 1 and 2 in R 1,2 (n) are the down-stream and up-stream sensor planes, N is the sample length, and f i (m) is the mth tomograms from ith planes of sensors. (5) is the general form of cross-correlation method for process tomography to seek for the correlation [27,28]. (5) is further adapted in [29] to become suitable for online calculation by updating R 12 (n) with the new kth tomograms, given by (6): where superscript (k) indicates that the related value is the kth sample. The major difference between (5) and (6) is that the former uses the whole input tomograms to produce mean pixel velocity distribution, whereas the latter does not, resulting in the latter presents better results. Therefore, the local velocity in this study is approximated by (6).
It is worth noting that due to the principle behind the crosscorrelation method and the nature of some flow regimes, i.e. stratified (wavy) flow and annular flow in horizontal pipeline or annular flow in vertical pipeline, it is extremely challenging to extract local velocity of the gas phase with the cross-correlation method. Therefore, the velocity for these flow regimes is estimated by the superficial velocity of gas phase in this paper. Figure 5 demonstrates an example of the calculation of IC resolution. Suppose an ERT system produces a set of electrical concentration tomograms with 20 × 20 mesh at the speed 1000 fps, pipe diameter is 100 mm, the minimum size This means that at temporal dimension, for the given conditions, during the time for the flow moving forward 20 mm, i.e. the maximal diameter of the small bubbles, 20 frames are collected by the ERT system. Therefore, every 20 tomograms need to be averaged as one frame in the temporally remeshed grid (the lower part of figure 5). It is worth pointing out that if the temporal dimension is less than 1, the actual distance the flow travels is larger than 20 mm for every two consecutive frames, and therefore interpolation is required to locate appropriate position for the flow floating 20 mm.
After the calculation of IC dimension and remeshing process, the original input tomograms are transformed to ICs (the bottom part of figure 5), and each IC physically reflects a 20 mm × 20 mm × 20 mm cube, within which a spherical bubble with maximal diameter 20 mm can be located in the next step.

Bubble generation
After the re-meshing, every IC is going to be processed to identify bubbles using a new lookup table (for small bubbles) and an enhance isosurface extraction algorithm (for large bubbles).

The new lookup table.
Unlike the conventional colourbased lookup table which maps scalar value to colours, the new lookup table is to transfer scalar values to small spherical bubbles. It is usually believed that measurement error is up to 5% owing to noise [8]. In other words, when the concentration is below 5%, it is assumed only conductive phase in the pipeline, which in our case is the liquid phase. In other words, the smallest size of a bubble is based on 5% concentration of the bubble in an interrogation cell. According to [30], the threshold value of the air concentration for the transition from bubbly flow to turbulent flow in a static bubble column with a diameter of 5cm was around 40% [31]. Thereafter, as an approximated value in the paper, it is assumed when volume fraction of gas is above 40%, large bubbles appear. Consequently, 40% is used as the threshold value beyond which an IC is assumed to be full of gas. When the concentration is between 5% and 40%, there would be a small bubble inside IC, of which the size is determined by the concentration. It is worthwhile to note that this study does not intend to improve the spatial resolution of electrical tomography. Instead, it tries to represent the concentration distribution of the gas phase from the reconstructed electrical concentration tomograms in bubble terms according to the above discussion.
For gas-liquid flow, concentration of an IC is calculated as: where c IC i is the gas mean concentration of the cell, v b j is the volume of ith bubble inside the given volume v IC i , and m is the number of bubbles inside the volume v IC i . In the context of gas-liquid flow, it is acceptable to approximate that when an IC is reasonably small, there is only one bubble inside an IC, i.e. m is 1, thereby (7) is rewritten to: where v b i is the volume of the bubble in ith cell. (8) can be further rearranged to: According to the dimension definition of IC, c IC i can be calculated by: By considering (3), (4), (9), and (10) altogether, the radius of the bubble with the volume v b i can be derived: Incorporating the assumptions into (10) as the transfer function to locate the entry of the new lookup table, and (11) as the calculation of bubble radius within an IC, the new lookup table can be established, as depicted in figure 6, in which T g and T l are the critical threshold values.
Further to the example in section 2.2, i.e. an IC physically reflects a 20 mm × 20 mm × 20 mm cube, the smallest size of a bubble in the cell can be derived accordingly to (11):

Small bubble mapping.
Once the bubble-based lookup table is built, it is ready to work on the interrogated concentration data. Figure 7 demonstrates an example of the application of conventional colour mapping and new bubble mapping to 2D stacked concentration distribution. In comparison to colour-based illustration, bubble-based one is able to reveal the size and position of small bubbles. As far as the large bubble in figure 7 is concerned, the outline is already drawn, but inconsistent with the real-world situation in terms of bubble shape. This is because no consideration of merging neighbour bubbles beyond particular size has been taken. In addition, figure 7 also illustrates the existence of some unreal bubbles in the red ellipse. In consequence, the next step is to merge neighbour bubbles to form larger bubbles.

Large bubble identification.
After the mapping, some cells may be filled with small bubbles, instead of the boundary segments during the interrogation. Those neighbouring cells with high gas concentration should be merged together to form a large bubble. After a cluster of the cells is located to form a large bubble, isofurface algorithm, e.g. Marching cubes [32], is applied to extract the interface between the bubble and liquid. Suppose a cell with the values of all vertices are larger than a given threshold, but this cell is located at a physical boundary, e.g. pipe wall, the algorithm, unfortunately, will not respect that there should have an isosurface segment as the physical boundary. The left image in figure 8 shows an example of an open contour, with a black contour of 0.4. It clearly demonstrates that without the consideration of the physical boundary, the contour becomes open, which is obviously inconsistent with reality. As a result, the algorithm has to be enhanced to reflect this situation. The right image in figure 8 illustrates an example of the application of enhanced isosurface algorithm with the consideration of the physical boundary.

Evaluation
The proposed approach is evaluated using concentration tomograms by ERT systems, on an industrial gas-liquid flow in both vertical and horizontal pipelines, covering almost all typical flow regimes in gas-liquid flow 3 . For each flow regime, three images are presented by high-speed camera, conventional colour mapping, and the proposed approach. All colour mappings of concentration are generated using the same colour mapping (from blue to red) and the same scale (0%, 100%). The local velocity of the gas phase was calculated by the aforementioned online cross-correlation method or the superficial velocity of gas phase when the method is inapplicable. When the method is applied, the cross-correlation length is fixed at 4000 frames, and the moving window is fixed at 500 frames. In order for the comparison to be more reasonable, all images for a given pipe orientation were generated to represent the same length, i.e. the distance of the flow passed through the pipe. Therefore, the total frames needed for generating the images reflecting the same travel distance of each flow condition can be derived by: where Dist is the travelling distance of the flow, DAS is the data acquisition speed of the ERT, and v i is the mean pixel velocity of the gas phase by the online cross-correlation method or the superficial velocity of the gas phase.

Gas-liquid flow in horizontal pipeline
A set of experiments was conducted in a horizontal pipeline with diameter of 200 mm on the flow testing facilities at TUV NEL, and the data was gauged by a commercialised ERT system, namely V5R [6], at the speed of 312.5 dual frames per Observed flow regime second (dfps). The mesh resolution of the 2D tomograms by V5R is 20 × 20. The actual travelling distance reflected by the images is 1.4 m. The selected flow conditions for the evaluation are listed in table 1, where flow reference of water superficial velocity, gas superficial velocity, and GVF (Gas volume fraction) were provided by the testing facilities. According to (12) and the actual flowing distance (1.4 m), the total frames for generating the images are listed in table 1. Since the local velocity of stratified (wavy) and annular flow is not derivable using the cross-correlation method, the superficial velocity is therefore applied to approximate the local velocity of the gas phase for these two flow regimes. Another point is that the error by measuring pure liquid is less than 3% during the experiment. Consequently, 3%, instead of 5% is chosen for T l . Given the above parameters/thresholds, visualisation results based on the bubble mapping are depicted in figure 9, along with the images from a high-speed camera taken through a photo chamber and the visualisation by conventional colour mapping. For a stratified flow in figure 9(a), the flow regime can be clearly recognised by all images. The situation changes when coming to a bubbly flow as shown in figure 9(b). The images by the high-speed camera and the colour mapping are incapable of visualising the small bubbles in the bubbly flow, whereas the bubble mapping shows how the small bubbles distribute in the pipe. For a plug flow (figure 9(c)) and a slug flow ( figure 9(d)), the camera and the colour mapping present similar results, i.e. the large bubbles can be located but the small bubbles are missing, whereas the bubble mapping complement the limitations by the methods. As far as annular flow is concerned, both the colour mapping and the bubble mapping are unable to demonstrate the thin water film at the top of the pipe. This is because the film is too thin to be identified by ERT system.

Gas-liquid flow in vertical pipeline
Another set of experiments was performed for an upward gasliquid flow in a 50 mm vertical section of the flow loop in the University of Leeds. 20 × 20 cross-sectional concentration data was obtained by another commercialised ERT system, called FICA [5], at the speed of 1000 dfps, and later stacked sequentially to form a 3D concentration dataset. The selected conditions for the evaluation are listed in table 2. In the images, the flow travels 0.4 m, and hence the total frames for each condition to generate the flowing distance are also listed in the table. As the aforementioned reason, the local velocity of the annular flow is replaced using the superficial velocity of the gas phase. 40% is used for T g and 5% for T l . The visualisation results, including bubbly flow, slug flow, and annular flow are presented in figure 10, accompanied with concatenated images by high-speed camera and colour-based rendering.  The resultant visualisation of a bubbly flow in a vertical pipe is different from the one in a horizontal pipe, as depicted in figures 9(b) and 10(a). In each group of images, e.g. figure 9(b), the top shows the photo from camera, the middle is the axial cross section of stacked tomograms, and the bottom shows the image produced by the bubble mapping approach. It is clearly demonstrated that both the camera and the bubble mapping can visualise the small bubbles, whereas the colour mapping does not. For a slug flow ( figure 10(b)), the flow regime can be recognized, but it is a little challenging to locate the exact size and distribution of the large bubbles due to the existence of small bubbles using the camera, whereas the stacked tomogram illustrates the approximate size and position of the large bubbles. In contrast, the bubble mapping reveals the small and the large bubbles. When it comes to an annular flow (figure 10(c)), both the camera and the tomogram can show the flow regime, but struggle to reflect the thickness of the liquid film surrounding the gas, whereas the bubble mapping can.

Impact of the critical parameters
Since the parameters, i.e. T g , T l , and IC dimensions, are critical for the approach, the impact of the values on the results is evaluated. For each parameter, three different values are selected, and corresponding visualisation results are presented, while    keeping the other parameters unchanged. Plug flow in horizontal pipeline (figure 9(c)) is chosen for the demonstration because it contains both small and large bubbles. Figure 11 depicts the influence of the threshold for pure liquid T l on the visualisation. The original rendering with T l = 3%, i.e. the top one in the figure, is set as the reference, and two more renderings with T l = 2% (the middle one) and T l = 8% (the bottom one) are generated. In the renderings, 174, 245, and 79 small bubbles are reconstructed, respectively, whereas the number of large bubbles is unchanged. In consequence, the threshold T l affects the number of reconstructed small bubbles, but does not affect the large bubbles, mainly because the threshold determines the proportion of liquid phase. That is, the higher the threshold is, the less the liquid is, and in turn the more the gas is.
The effect of the threshold T g is presented in figure 12, in which three different values, including 40%, 50%, and 60%, are applied, and the outcomes are from top to bottom respectively. From visualisation point of view, the change does not have a distinguishable effect on the outcomes. The number of acknowledged small bubbles also proves it, in which there are 174, 174, and 183 small bubbles, respectively. Since T g determines the full occupation of an IC by gas, it is crucial to the volume of large bubbles, and the boundaries between large bubbles and liquid, although the figure does not explicitly illustrate it.
Another important parameter is the dimensional definition of IC, and figure 13 shows the comparison. In essence, the size of IC is vital for the size of small bubbles, since the volume is not beyond an IC. In addition, it also influences the number of mapped small bubbles, since smaller IC results in more amount of IC with a given concentration distribution. Consequently, the smallest IC results in the bubbles of the smallest in size and the most in number, as shown by the middle one in figure 13, in which totally 3251 small bubbles are located with the IC of 2 × 2 spatial dimension, compared to the ones by the IC of 4 × 4 and 5 × 5 with 174 and 158 small bubbles, respectively. In contrast, the appearance of large bubbles is irrelevant to the parameter.

Quantitative analysis of plug flow
Four cross-sectional images at different positions at the plug flow (figure 9(c)) are chosen for quantitative comparison. They are a large bubble, a large bubble with small bubbles, tail of a large bubble with small bubbles, and small bubbles. The images are demonstrated in figure 14, and mean air concentrations derived from the images are presented in table 3. It is worth noting that the cross-sectional images from the bubble mapping result includes pipe wall, but the corresponding concentration is derived excluding the wall. When there is only one large bubble, the images by ERT and bubble mapping (figures 14(a) and (e)) are similar, and the mean gas concentrations are very close (Case 1 in table 3). This is mainly because of the application of SBP reconstruction which makes the boundary between the large bubble and liquid sharp. When small bubbles are introduced, both images and mean gas concentration are different. As far as mean gas concentration is  concerned, bubble mapping results tend to underestimate it. This is primarily because of the assumption of noise, i.e. mean concentration less than 5% is assumed due to noise and therefore ignored.

Conclusion
A novel approach, namely bubble mapping, has been proposed to visualise gas-liquid flow based on electrical concentration tomograms. With the reasonable assumptions, small bubbles can be reconstructed and large bubbles can be located by a bubble lookup table. Small bubbles are randomly located in a cubic volume, and their sizes are fractional to the mean concentrations of air in the cubes. In contrast, large bubbles are random in shape and location, reconstructed by a cluster of neighbouring cubes with the local mean concentrations beyond a threshold value and an enhanced isosurface approach. The evaluational results demonstrated that through this approach, colour-based flow phase distributions are transformed to bubble-based visualisation, which is capable of significantly improving flow regime visualisation and visual recognition, and thus eliminating the visual ambiguities caused by the conventional colour mapping-based approaches. In addition, the bubble mapping results are able to convey more detailed information which the colour mapping ones are not, e.g. the bubble size increases when they are closer to the upper wall of the pipe in horizontal flow (e.g. figure 9(b)). It is worth noting that when large bubbles are involved, the overall mean concentration is not conserved any more due to the threshold value selected is 0.4, rather than 1.0, which may result in the overestimation of large bubble size. However, it has little influence on the visualisation and visual recognition of flow regimes. Despite of the feasibility of the proposed approach, there are still a few aspects to be addressed in the future studies. The threshold value used for classifying small bubbles and large bubbles in the bubble lookup table needs to be studied for the theoretical support, e.g. determining the value with an inverse method [33]. Since the visual results are dependent on the quality of input data, i.e. the spatial resolution of electrical concentration tomograms, the resultant visualisation could be improved by applying advanced algorithms in both tomography and computer graphics, such as the one by [15,34], respectively. However, the challenge still remains in overcoming the underestimation error from imaging of large non-conductive object with the current status of ERT. It is also worthwhile to indicate the potential artefact introduced, particularly the small bubbles appeared around a distinctive large bubble, which is caused by the point-spread effect of the tomographic imaging although they are not playing a predominant role in flow regime visualisation.