Multivariate Image Analysis (MIA) for Industrial Flare Combustion

A new approach for flare monitoring is proposed so that flare combustion efficiency can be predicted online in industrial plants. Multivariate image a...
0 downloads 0 Views 3MB Size
Article pubs.acs.org/IECR

Multivariate Image Analysis (MIA) for Industrial Flare Combustion Control David Castiñeira,† Blake C. Rawlings,‡ and Thomas F. Edgar* Department of Chemical Engineering, University of Texas at Austin, 1 University Station C0400, Austin, Texas 78712-023, United States ABSTRACT: A new approach for flare monitoring is proposed so that flare combustion efficiency can be predicted online in industrial plants. Multivariate image analysis (MIA), which is based on principal component analysis (PCA) and projection to latent structures (PLS), has been applied to flare combustion systems in order to predict their resulting combustion efficiencies, as a function of the crosswind velocity, using simulated results, and as a function of steam or air flow rates, using experimental tests of a full-size flare. The results show that a multivariate regression model based on flare color images can be used to predict the flare performance over a range of operating conditions for steam-assisted flares. Therefore, simple two-dimensional color images of industrial flares may be a fast, accurate, and inexpensive approach for online monitoring of these industrial combustion systems. This would allow for developing effective flare control and mitigation strategies. that the flame properties and resulting combustion efficiency can be accurately predicted online from flare images. As a starting point, we propose the implementation of MIA to predict the resulting combustion efficiency of industrial flares subject to different crosswind velocities and different steam or air injection rates. Note that most chemical and petrochemical plants typically use an assumed 98% combustion efficiency in calculating emissions from flares, regardless of the amount of steam or air injection and crosswind velocity. However, several studies have shown that flares subject to high crosswinds may operate at low efficiencies, resulting in significant hydrocarbon emissions to the atmosphere.8−911 This work shows how moreaccurate combustion efficiencies of flares subject to different operating conditions (including crosswind speeds) can be estimated in real time from flare image analysis.

1. INTRODUCTION Industrial flares are units designed to safely dispose of the waste gases from chemical and petrochemical plants. In these units, a waste gas stream is fed through a stack anywhere from 10 m to over 100 m tall and is combusted at the tip of the stack, leading to a large turbulent flame. Flares have been shown to be a potential source of hydrocarbon gas emissions to the atmosphere, so there is a recent need to develop control and mitigation strategies to minimize their impact on air quality. Any flare control strategy would require some type of measurement of the flare performance. Unfortunately, for industrial flares, it is very difficult to obtain online efficiency and emission measurements for several reasons. First of all, exact estimation of the flare efficiency would require collecting the entire plume of combustion gases, which is something that is impractical in most industrial plants. On the other hand, taking local measurements within the flame is complicated by the large size of the flames, their high turbulence, and the variable ambient conditions (crosswind, rain, etc.) that exist in these industrial units. Under these circumstances, an inferential sensor based on flare image analysis seems very attractive for practical implementation. The main advantage of using image analysis to monitor a combustion process is based on the fact that different colors within the flame indicate different local temperatures and, consequently, different combustion regions. Thus, multivariate statistical analysis of the flare image colors could be a valid approach to predict the overall flare performance. Principal component analysis (PCA) and partial least-squares (PLS) can be applied for this purpose, so most of the mathematical analysis is restricted to the latent variable feature space, rather than in the observed image space. Different visualization methods have been proposed over the last decades for combustion systems.1−7 However, none of these studies has focused on the application of multivariate image analysis (MIA) to industrial flares. Here, we investigate the applicability of MIA for flare monitoring and control, so © 2012 American Chemical Society

2. MULTIVARIATE IMAGE ANALYSIS (MIA) FUNDAMENTALS Color images represent a class of multivariate image.12 A multivariate image is a stack of “congruent” images, that is, images where pixel positions match for all images in the stack. Multivariate imaging has been introduced in different scientific fields, where congruent images of the same object or scene have been traditionally obtained by using different radiation energies or wavelengths. Multivariate images can also be obtained by combining images resulting from different instruments or even by combining a univariate image with copies of it that are derived from filtering operations. An RGB color image is a multivariate image with three variables (e.g., the horizontal dimension, the vertical dimension, and a third dimension, which is given by the R, G, and B channels). Special Issue: Industrial Flares Received: Revised: Accepted: Published: 12642

February April 18, April 18, April 18,

7, 2012 2012 2012 2012

dx.doi.org/10.1021/ie3003039 | Ind. Eng. Chem. Res. 2012, 51, 12642−12652

Industrial & Engineering Chemistry Research

Article

regression model that directly relates the flare efficiency with the flare image. The partial least squares (PLS) regression technique models the relationships between observed variables by means of latent (not directly measured) variables. It constructs predictive models to analyze large multivariable collinear data sets. The method aims to identify the underlying factors, or linear combination of the X variables, which best model the Y dependent variables. Figure 1 gives a schematic outline of the method.15

A typical multivariate image usually contains an enormous amount of data (for a single 512 × 512 × 3 (RGB) color image, a total of 786 432 pixel intensities are needed to describe the image information). However, these data are usually highly correlated, because of the overall brightness or darkness ranges throughout the image, which overrides the color differences, so pixel variations within the raster are similar, despite the band. Principal component analysis (PCA) has been utilized by Yu and MacGregor6 to compress such data and project it onto a reduced dimensional subspace through linear combinations of the original multivariate data. This step is explained in more detail in the next section. For the image plane (I × J dimension), the planar geometry can be converted to a single dimension N, where N = I × J, and K is on the three-dimensional RGB intensities for a given pixel. ⎡ R1,1 ⎢ ⎢ ⎢ ⎢ R1, J X=⎢ ⎢ R 2,1 ⎢ ⎢ ⎢⎣ RI , J

G1,1 B1,1 ⎤ ⎥ ⋮ ⎥ ⎥ G1, J B1, J ⎥ G2,1 B2,1⎥⎥ ⎥ ⋮ ⎥ GI , J BI , J ⎥⎦

Figure 1. Schematic representation of the indirect modeling for partial least squares (PLS).

PLS creates orthogonal score vectors by maximizing the covariance between different sets of variables.16 The predictor and predicted (response) variables are each considered as a block of variables. Consider ? ⊂ N to be an N-dimensional space of variables representing the first block; similarly, consider @ ⊂ M to be a space representing the second block of variables. After observing n data samples from each block of variables, PLS decomposes the (n × N) matrix of zero-mean variables X and the (n × M) matrix of zero-mean variables Y into the form:

Given a two-dimensional image matrix X of dimension (N × K), PCA decomposes this matrix into simpler matrices Ma: A

X=

∑ Ma

(1)

a=1

The smallest value of A for which eq 1 is correct is called the rank of X. The Ma matrices all have dimensions of (N × K). The Ma matrices also have the special property that all have rank = 1. Therefore, they can be represented as the outer product of two vectors, t and p:

Y = UQT + F X = TPT + E

where T and U are (n × p) matrices of the p extracted score vectors; and the (N × p) matrix P and the (M × p) matrix Q represent matrices of loadings; and, the (n × N) matrix E and the matrix F are the matrices of residuals. The difference between partial least squares (PLS) and principal component analysis (PCA) is that PLS does not only try to minimize the residuals E and F, but it maximizes the correlations between U and T. The PLS method, based on the nonlinear iterative partial least squares (NIPALS) algorithm,17 finds the weight vectors w and c such that

A

X=

∑ tapa a

(3)

T

(2)

where A represents the number of principal components, ta represents the score vectors of size (N × 1) (these vectors are orthogonal), and pa represents the loading vectors of size (K × 1) (these vectors are orthonormal). Because an RGB color image is, by definition, a three-way array (of dimensions I × J × K), it is convenient to first transform this array into a two-way array or matrix (dimensions N × K) so then PCA can be applied. For RGB color images, the maximum number of principal components is 3, so the PCA decomposition will result in three ta vectors and three pa vectors, one for each component. The first component explains the greatest amount of variance in X, the second component explains the next greatest variance, and so forth. The number of principal components necessary to extract most of the meaningful information can be determined by several methods.13,14 Typically, the first two components explain most of the variability of the original data in color images. Therefore, PCA represents a very useful tool for data reduction and data interpretation in the analysis of flare color images. Furthermore, different image features can be easily extracted from the resulting score plots that have more physical meaning. Those image features then can be regressed against some known values of the flare efficiency. This would create a

[cov(t , u)]2 = [cov(Xw , Yc)]2 = max |r|=|s|= 1[cov(Xr , Ys)]2 (4) T

where cov(t,u) = t u/n denotes the sample covariance between the score vectors t and u.

3. WIND TUNNEL SIMULATION STUDY This work analyzes the application of MIA to a series of natural gas flares located in a closed-loop wind tunnel. The actual experimental flares were studied at the University of Alberta and the National Research Council of Canada. Measurement of experimental combustion efficiencies at different crosswind velocities were reported by Johnson and Kostiuk.10 For these flares, the jet exit velocity of the fuel (Vj) was held approximately constant at 2 m/s and the crosswind speed (U) was varied from 1 m/s to 11 m/s. The external cold-flow Reynolds number (Re) ranged from 1570 to 17 270 as the 12643

dx.doi.org/10.1021/ie3003039 | Ind. Eng. Chem. Res. 2012, 51, 12642−12652

Industrial & Engineering Chemistry Research

Article

crosswind increased from 1 m/s to 11 m/s. Only limited experimental data were available for the flares (e.g., image data), so it was decided to simulate these wind-tunnel flares using computational fluid dynamics (CFD),18 in order to generate representative images. In the simulation work, a set of threedimensional (3-D) computational models were developed and suitable models were applied to simulate the complex combustion phenomena and flame downwash. Large eddy simulation (LES) was implemented to model turbulence, so the turbulent fluctuations of the buoyant plume and the unsteady flame shape dynamics could be captured. In fact, LES resolved the large length and time scales responsible for controlling the dynamics of the flame, with models for more homogeneous smaller scales. Unfortunately, resolving even the large length and time scales requires large amounts of computational resources. Therefore, this task had to be done in a massively parallel environment with the aid of computational tool that was scalable. The LES simulations were carried out at the Combustion and Reaction Simulation (CRSIM) research group at the University of Utah, which has developed a massively parallel LES code (ARCHES). This code is linearly scalable up to 1000 processors. A total of six flare cases were analyzed in this simulation work: cases B, C, D, E, F, and G, as named by the authors of the actual flare lab experiments.10 The simulation results for the flame temperature were displayed in a two-dimensional (2-D) plane using a volume rendering approach, which allows the visualization of 3-D data by computing 2-D projections of a colored semitransparent volume. Results showed that flare images obtained from LES simulation were in very good agreement with the actual photographs of the experimental flares. One advantage of using simulation images for MIA is that they do not contain noise or complex backgrounds that may make the analysis of the flame difficult. In fact, the computational images used in this work did not require any special preprocessing. For each of the six flares considered in this study, 20 color images were retrieved that were 0.1 s apart. Hence, a total of 120 images were analyzed. Because of the inherent turbulent fluctuations of the flame, flare images changed significantly in very short time, even for a constant crosswind velocity. Figure 2 shows three sample images for each of the flare cases analyzed here, where the images are given at a reference time (t0), at one second (1 s) apart (t1), and at two seconds (2 s) apart (t2). Notice that the flare stack has been conveniently removed from the images to facilitate the subsequent image analysis.

Figure 2. Sample images of the flares taken at three consecutive seconds (t0, t1, and t2). Images are shown for increasing crosswind velocities (cases B, C, D, E, F, and G).

in the image to remove all these background features, as discussed later in this paper. (2) Unfold the three-way arrays (X) into two-way arrays or matrices (X). there are two practical ways to unfold a 3-d array into a 2-d array or matrix.12 the resulting matrices are long matrices (number of rows = number of pixels) with only three columns. (3) Apply principal component analysis (PCA) to decompose X into loading and score vectors. A multivariate image typically contains an enormous amount of data that is highly correlated. PCA can be used to project it onto a reduced dimensional subspace through a few linear combinations of the original multivariate data. The vector of values of the principal component at each pixel location in the scene space is given by the score vector ta, where

4. PROCEDURE FOR IMPLEMENTING MIA The practical implementation of multivariate image analysis (MIA) followed a procedure that is similar to the approach of Yu and MacGregor,6 which was originally developed for analysis of industrial boilers. The procedure allows creation of a regression model that can be used to predict the flare combustion efficiency at a given crosswind velocity. (1) Convert color images into three-way arrays (X) of data according to RGB intensities. Note that the analysis of real flare images would be more complex due to the presence of the flare stack and a potential interfering background (for example, a highly luminous sky may affect the color analysis of the flare image). Under these circumstances, some sort of processing and transformation operations would be performed

ta = Xpa

(5)

To perform PCA decomposition of a matrix, a kernel algorithm19 can be used where SVD is performed on the kernel matrix XTX to obtain the loading vectors of the original matrix directly. For color images, the dimensions of XTX will be (3 × 3), so the computational work is greatly simplified. Furthermore, since we are dealing with a set of images, the kernel matrix must be calculated as the summation of the covariance matrices of the individual images, that is, as ∑jXjTXj. Once the loading vectors of the kernel matrix have been obtained, the corresponding score vectors for the individual 12644

dx.doi.org/10.1021/ie3003039 | Ind. Eng. Chem. Res. 2012, 51, 12642−12652

Industrial & Engineering Chemistry Research

Article

t1 and t2 values; low-temperature points within the flame would also produce similar t1 and t2 combinations. On the other hand, extinguished points or background points would produce significantly different t1 and t2 values. Consequently, we can plot the t1 and t2 score vectors against each other (scatter plot), so that all the points in the image with the same spectral characteristics would plot on top of each other or at least in the same neighborhood in this score plot. Hence, score plots are very useful in the analysis and monitoring of online flares images, since they simplify the extraction of combustion (reactive) regions from the images. Since an image typically has a huge amount of pixels, many t1 and t2 values may overlap in the scatter plot, making it difficult to identify those points of interest. A common approach in multivariate image analysis is to define a density score plot, where an intensity value is assigned to each point in the s1−s2 domain. The intensity assigned to each point (s1, s2) is equal to the number of pixels in the original image that yielded the scores s1 and s2. The densities can be saved in a 256 × 256 density matrix where each entry is a density value. This density matrix can be easily plotted as an image in order to identify clusters of points as shown in Figure 4. This allows one to identify those points associated to the combustion region, background, noise, etc. Mathematically, the density plot is given by

images can be calculated from eq 2, leading to the singular values and percentage of variability shown in Table 1. Table 1. Singular Values of the Kernel Matrix and Percentage of Variability Explained by Each Principal Component 1st component 2nd component 3rd component

λ1

percentage of variability explained

7.438 × 1010 0.336 × 1010 0.025 × 1010

95.38% 4.31% 0.32%

Because the two first components explained 99.69% of the variability of the original data, the third component can be removed without losing significant information. In Figure 3, a

Figure 3. Reconstructed image from the first two components (left side) vs original flare image (right side) of Flare B.

TTij =

∀ p , slp = 1, s2p = j

i , j = 0 , ..., 255

p

reconstructed image (left side) from these first two components is compared to one of original images of Flare B (right side). This image is basically a color representation of the matrix of intensities obtained through eq 6:

(7)

The density score plot in Figure 4 displays higher density values with a brighter color. The colors scale so that a point (s1,s2) that does not match any pixels in the original image (a density of zero) is black, and the point with the highest density is white. This plot clearly differentiates two large clusters of points (the left arc and the right arc), which correspond to the black background and the flame region of the image respectively. In the next section, a masking procedure is described that helps in extracting those features of interest from these density score plots. It is worth noting that flare images at the same combustion condition produce similar density score plots. This is clearly shown in Figure 5, where two density score plots (corresponding to images obtained 2 s apart) are shown for flare cases B, C, D, E, F, and G, respectively. As we can see from this figure, the resulting density score plot is very similar for flare images obtained under the same crosswind conditions,

[Reconstructed image]ij = [t1p1T + t 2p2 T ] i = 1..., I , j = 1, ...K

∑1

(6)

As can be seen from this figure, the first two components summarize the dominant spectral features of the image; therefore, the original image is almost entirely reconstructed from these components. From this point onward, a much simpler analysis of the image can be done in terms of the first and second scores (t1 and t2, respectively). (4) Obtain the score plots and density score plots. Regardless of their spatial location, similar occurrences or features in an individual image would produce the same combination of score values (t1 and t2). For example, all the high temperature points within the flame would produce similar

Figure 4. Score plot (left side) and density score plot (right side) of the s1 and s2 score vectors. 12645

dx.doi.org/10.1021/ie3003039 | Ind. Eng. Chem. Res. 2012, 51, 12642−12652

Industrial & Engineering Chemistry Research

Article

Figure 5. Density score plot of s1−s2 pair for flare cases B, C, D, E, F, and G. 12646

dx.doi.org/10.1021/ie3003039 | Ind. Eng. Chem. Res. 2012, 51, 12642−12652

Industrial & Engineering Chemistry Research

Article

Figure 6. (a) Mask region over the density score plot; (b) original image for Flame B; and (c) flame luminous region decided by the mask (all points outside the mask are set to a gray color).

Figure 7. (a) Mask region over the density score plot; (b) real flare image; and (c) flame luminous region isolated by the mask (all points outside the mask are set to a gray color).

regardless of the turbulent and unsteady fluctuations of the buoyant plume. On the other hand, for flare images obtained at different crosswind velocities, the locations of pixels in the score plots changed noticeably. The differences would be more noticeable in the case of real flare images. Thus, the flare combustion efficiency can be estimated based only on fundamental analysis of the dominant spectral features of the image and regardless of the turbulent fluctuations that may produce significant changes in the flame shape and image. This is the key result that makes multivariate image analysis of flare images a useful tool to monitor the combustion efficiency. (5) Extract features from score plots. The goal of studying the density score plots was to be able to identify those clusters that contained meaningful pixel classes. Typically, a few clusters are observed in a given image, where one cluster may represent, for example, the image background while another cluster may represent the high-temperature regions of the image (e.g., the actual flame). Furthermore, in many images, there are pixels that do not belong to any well-defined cluster, but have fuzzier cluster membership. In this situation, it becomes important to correctly delineate the pixels of interest by applying multivariate segmentation. Several multivariate segmentation techniques have been proposed by Esbensen and Geladi20 and Geladi and Grahn.12 The basic idea is to extract the flame luminous region by choosing a mask in the score plot space. The boundary of the mask can be easily obtained by a trial-and-error process of drawing a polygon and isolating the pixels inside that polygon. We then can highlight the associated pixels in the image space, and iterate until one obtains a mask that extracts the feature of interest. For this particular work, the masking procedure was rather simple, because of the nature of the images studied (e.g., images were obtained directly from computer simulation, so the

different clusters could be easily identified). For practical implementation in industrial flares, a novel automatic masking procedure could be applied21 that helps delineate the pixels of interest. Figure 6 shows the ability of a mask to extract the flame’s luminous region for a sample image of Flare B. The mask drawn over the density score plot is shown in Figure 6a. If all pixels falling outside this mask in the image space (e.g., in the original image, Figure 6b) are set to have a gray color, then Figure 6c is obtained. Clearly, the mask is capable to isolate those pixels corresponding to the flame’s luminous region. For industrial flare images, the masking procedure would be more difficult, because of the presence of a luminous sky in the background or other objects around the flares. However, the flame’s luminous region can still be isolated by choosing the right mask in the score domain, as shown in Figure 7 for a sample image of a real flare. As we can see from this figure, an effective flare MIA can also be performed for complex backgrounds that may exist in the industrial plant. Once pixels corresponding to the flame’s luminous region have been delineated, the user is in position to extract useful information from the image. In fact, certain features or “image measurements” can be easily extracted from the masked score plot that have physical meaning. Those measurements will then be regressed against the experimental combustion efficiencies of the flares, so that a predictive model is obtained. First, a binary matrix M that fits over the 256 × 256 score plot as a matrix, which has an element equal to 1 if the element location lies under the luminous flame mask and 0 if it does not, must be defined. Hence, the following feature variables will be used to build the regression model: Area of the Luminous Region (A). The area of the flame image that corresponds to the luminous region can be easily 12647

dx.doi.org/10.1021/ie3003039 | Ind. Eng. Chem. Res. 2012, 51, 12642−12652

Industrial & Engineering Chemistry Research

Article

Figure 8. Flame luminous region (A) and average location of the s1 pixel in the score plot (s1m) for the different flare cases analyzed in this work. These feature variables are obtained from the 120 flare images considered here.

the flame is captured in the image. Therefore, the brightness of the nonluminous area presumably gives information on the length and/or the volume of a flame. We compute the average brightness of the nonluminous area by

computed by counting the number of pixels falling into the flame mask. Mathematically, A=

∑ TTij

∀ (i , j), Mij = 1 (8)

i,j

where M is the binary mask defined in the previous section. Flame Brightness (F). Flame brightness can be obtained by integrating the luminous intensity level contributed from all pixels falling inside the luminous region. The luminous intensity level for any location in the score plot is computed by converting the color plane obtained by eq 10 to a gray-scale plane with elements corresponding to the luminous intensity at each location. Each location in the s1−s2 score plot represents a certain type of color (ignoring the residual variation in t3) [ R G B ]ij = t1(i)p1T + t 2(j)p2 T

W=

s1m =

t 2(j) =

255 255

+ t1min

⎡ 0.299 ⎤ ⎢ ⎥ Lij = [ R G B ]ij ⎢ 0.587 ⎥ ⎢⎣ 0.114 ⎥⎦

Nc = (10)

N

Uniformity of Flame Brightness (U). The uniformity of the flame brightness is defined as the standard deviation of the flame brightness throughout the luminous region, ∑i , j TTijL2 i , j − ∑i , j TTijLi , j

i , j = 0 , ..., 255

A

s 2f =

∑i , j TTijj A

∀ (i , j), Mij = 1

=

∑1

∀ (i , j), (TTijMij) ≠ 0 (16)

5. REGRESSION MODEL FOR SIMULATED FLARE DATA A total of 120 images were used for multivariate image analysis (MIA) of the simulated data. For each image, nine feature variables were obtained according to eqs 8−16. Two of these feature variables (flame luminous region (A) and the average location of s1 pixels in the score plot, s1m) are shown in Figure 8 for the different flare cases analyzed in this work. The trend followed by these two variables is exactly the opposite. As expected, increasing the crosswind velocity reduces the flame luminous region, and a similar trend is observed for other parameters, such as the flame brightness (F). The trend observed for the average s1 pixel location, s1m, and some other extracted variables is not so obvious, from a physical point of view, since they are more related to the spectral information of the images. Most importantly, all the parameters defined in this work seem to have a strong correlation with the crosswind velocity and, consequently, with the resulting combustion efficiency.

(11)

∀ (i , j), Mij = 1

∑i , j TTiji

i,j

∀ (i , j), Mij = 1

∑i , j TTij

∑i , j TTijj

Total Number of Colors Appearing in the Flame Luminous Region (Nc). The number of the different colors appearing in the flame region is the area of the flame region in the score plot space:

i , j = 0 , ..., 255

i,j

U=

N

s 2m =

(15)

where the last column gives the conversion coefficient vector that converts the element RGB signal into a corresponding luminous intensity. The overall flame brightness is then calculated as

∑ TTijLi ,j

∑i , j TTiji

+ t 2min

Given the color plane of an image defined by eq 9, the luminous intensity level for location (i,j) is given by matrix Lij, which is defined in eq 10:

F=

(13)

Average Color of the Flame Luminous Region (s1f, s2f). The average color of the flame’s luminous region is defined as the average location of the pixels belonging to the flame luminous region: s1f =

j(t 2,max − t 2,min)

∀ (i , j), Mij = 0

(14)

i , j = 0 , ..., 255

where i(t1,max − t1,min)

∑i , j TTij

Average Color of the Entire Flame Image (s1m, s2m). The average color of the entire flame image is expressed as the average location of pixels in the score plot:

(9)

t1(i) =

∑i , j TTijLi , j

∑i , j TTijL2 i , j − B A (12)

Average Brightness of the Nonluminous Area (W). The flame image is a 2-D projection of a 3-D field and only part of 12648

dx.doi.org/10.1021/ie3003039 | Ind. Eng. Chem. Res. 2012, 51, 12642−12652

Industrial & Engineering Chemistry Research

Article

not available currently. The regression models can be extended to predict efficiencies for steam-assisted and air-assisted flares as discussed in the next section.

Before proceeding with the partial least-squares (PLS) analysis, some form of image filtering should be performed to reduce the variability of the extracted data that is due to the highly turbulent nature of the flames. Thus, a time-average filter with an average moving window with length 4 was applied. Finally, PLS regression was performed using the nine feature variables obtained from the flare images as predictors. These variables were regressed against the experimental combustion efficiency observed for these flares.10 The NIPALS algorithm was applied. Half of the observations were used as training set, and the other half were used as a test set. Note that other training/testing schemes are explored in Section 6 with actual (as opposed to simulated) data. Figure 9 shows a comparison between the predicted and the observed combustion efficiencies.

6. FULL-SCALE FLARE TEST RESULTS In May 2009, the Texas Commission on Environmental Quality (TCEQ) funded the University of Texas at Austin (UT Austin) for a Comprehensive Flare Study project to conduct field tests to measure flare emissions and collect process and operational data in a semicontrolled environment to determine the relationship between flare design, operation, vent gas lower heating value (LHV) and flow rate, destruction and removal efficiency (DRE), and combustion efficiency (CE). The TCEQ’s primary objectives for this study were • to assess the potential impact of vent gas flow rate turndown on flare CE and volatile organic compound (VOC) DRE; • to assess the potential impact of steam/air assist on flare CE and VOC DRE under various operating conditions, including low vent gas flow rates; • to determine whether flares operating over the range of requirements stated in 40 Code of Federal Regulations (CFR) §60.18 achieve the assumed hydrocarbon DRE of 98% at varying vent gas flow rate turndown, assist ratios, and vent gas heat content; and • to identify and quantify the hydrocarbon species in flare plumes visualized with passive infrared cameras. To make the results of this study most directly applicable to industrial-scale operations, the field tests performed for this study were conducted on full-scale industrial design flares. Specifically, the flare designs selected were the John Zink Models EE-QSC-36′′ Flare Tip (36-in. diameter) with three pilots (air assist), with maximum capacities of 937 000 lb/h and 144 000 lb/h, respectively. These sizes and design configurations were selected so that they represented a large number of flare models currently in the field. To focus on low LHV (lower heating value) vent gas streams and still comply with 40 CFR §60.18, a LHV of 350 Btu/scf ±80 Btu/scf was selected as the lowest target LHV for the vent gases used in the test plan. To obtain additional data on the effect of LHV on DRE and CE, a second feed gas with LHV of 600 Btu/scf ±80 Btu/scf was also included in the test plan. The test plan consisted of a matrix of flare operating conditions designed to provide data that would be the basis to address as many of the study objectives as possible. This matrix of operating conditions included two low vent gas flow rates for

Figure 9. Comparison between the predicted and observed combustion efficiencies for increasing crosswind velocities. Predictions were obtained through partial least-squares (PLS) regression.

Figure 9 shows that the combustion efficiencies can be satisfactorily predicted from the regression model based on the feature variables extracted from the flare images. The observed experimental values were measured as average values for constant crosswind velocities, so they present a plane profile in Figure 10. Predicted values show some oscillation, which is partially due to the flame turbulence and the characteristic pulsating flow. In any case, good estimates of the combustion efficiency were observed. Clearly, these results indicate that the spectral information contained in the flare image is strongly related to the thermochemical state of the flare and the resulting combustion efficiency. This approach for estimating flare combustion efficiencies could be of great advantage for online monitoring of flares since quick estimates of the combustion efficiencies are

Figure 10. Combustion efficiency (CE) versus steam assist for all test series S3 and S4 for gas with LHV = 350 Btu/SCF. 12649

dx.doi.org/10.1021/ie3003039 | Ind. Eng. Chem. Res. 2012, 51, 12642−12652

Industrial & Engineering Chemistry Research

Article

the steam flare (937 and 2342 lb/h) and two low LHVs (300 and 600 Btu/scf). For the air-assisted flare, 359 and 937 lb/h vent gas flow rates and the same two low LHVs used for the steam flare were used. The vent gas composition used was a 1:4 ratio of Tulsa natural gas to propylene diluted to achieve the desired LHV. Air and steam assist rates used varied from the amount used to achieve the incipient smoke point to an amount near the snuff point. All of the tests in this study were conducted under conditions that are in compliance with all criteria of 40 CFR §60.18. All operating parameters for the flare were measured and monitored during each test run. The CE and DRE of the flare for each test point were determined by continuously extracting a sample from the flared gas beyond the point in the plume where all combustion had ceased and then analyzing the sample at a rate of 1 Hz, using a suite of analytical instruments operated by Aerodyne Research, Inc. A carbon balance was performed on the constituents in the sample, compared to the constituents in the vent gas flow, and the appropriate quantities were used to calculate DRE and CE. Two remote-sensing technologies were also employed in the study and have been compared to the extractive measurement results. Figure 10 illustrates the significant effect of steam addition to combustion efficiency. It is clear that as excess amounts of steam are added, relative to the feed gas rate, the combustion efficiency drops off precipitously to CE values well below 90%. There is an operating range where CE values over 98% can be achieved. However, because the flow rate or heating value of a flare gas can decrease, then oversteaming or increased air injection to the flare can occur if the flare operator does not make compensating changes in the rate of air or steam addition. Without appropriate measurements on flare gas flow or composition or an online estimate of CE, then there is no information upon which the operator can act. MIA was used to assess the relative impact on combustion efficiency by operating variables such as vent gas flow, steam or air assist, flame temperature, and the presence of certain VOCs. The resulting PLS model could also be used to better understand the performance data obtained in the flare tests and the effect of such parameters as crosswind, vent gas flow rate and composition, and air and steam assist at operating points that were not run in the tests. From an image of the flare, the set of feature variables extracted were based on the appearance of the flame. These variables included the area, angle of the flame (slope), brightness (intensity), and color (hue). The image background was easy to screen out. Assuming that the appearance of a flare is related to its CE, these values were then used to generate a model to predict the CE. The regression results are recorded both as the predicted CE from each input image, and also as the average predicted CE for the entire duration over which the flare was recorded. The rationale behind using the average CE value over the entire recording is that flares are not steady processes; instead, they continually change shape and appearance. Therefore, a single image cannot be used to accurately represent the state of the flare, and an average over a reasonable time-scale (on the order of minutes) must be used. The John Zink-UT tests that resulted in a visible flame are briefly described in Table 2. The other tests either did not produce visible flames (so image analysis in the visible spectrum was not possible), or did not differ appreciably in appearance from those flames described. The descriptions are provided to

Table 2. (a) Test Descriptions and (b) Combustion Efficiency and Operating Conditions of the Different Tests (a) Test Description test

appearance

Bright, small flame. Roughly 45° angle. Longer, dark orange flame. 45° angle. Slow movement in the wind. S4.1 Very faint flame. S5.6 Weak flame, almost horizontal. S6.1 Similar shape to S1.5, but not as bright. In some frames, parts of the flame are cut off by the camera. A2.1, A2.3, Bright, small flame. Flame does not extend much past the flare A2.4, tip, and is not long enough to be significantly bent over by the A2.5 wind Bigger, darker orange, less vigorous flame. Bent over nearly A3.1 horizontal in the wind. A3.6 Weaker flame than A3.1, otherwise similar. A6.1 Bright, small flame, bent nearly horizontal in the wind. A6.4, A6.5 Vary between almost invisible and a flame similar to A2.1 (b) Combustion Efficiency and Operating Conditions S1.5 S3.7

test

combustion efficiency, CE (%)

crosswind (mph)

wind directiona (deg)

lower heating value, LHV (Btu/scf)

assist rateb (lb/h)

S1.5 S3.7 S4.1 S5.6 S6.1 A2.1 A2.3 A2.4 A2.5 A3.1 A3.6 A6.1 A6.4 A6.5

99.9 99.2 95.0 92.6 99.2 95.5 94.0 87.6 91.8 99.1 88.9 99.3 86.3 81.6

8.0 7.1 5.6 9.6 8.8 12.8 10.1 10.0 13.3 10.3 11.9 15.9 14.1 15.5

176 133 146 180 180 174 180 185 174 178 170 180 180 185

2145 346 350 590 609 2125 2108 2113 2124 339 338 584 585 584

3794 228 536 463 1003 83818 88791 148799 119580 19387 60121 11404 40584 56594

a

A wind direction of 180° is sideways, from the perspective of the camera used to record the flares. bAssist rate refers to the flow rate of the assist gas in the steam-assisted and air-assisted flares.

indicate how much the appearance of the flare varied between tests, as image analysis relies on variations in the appearance of the flare in order to predict the combustion efficiency accurately. The operating conditions used in the regression also are shown in Table 2. The image analysis code was applied to the usable tests (those with a visible flame) for both the air-assisted and steamassisted flares. For each test, five minutes of image data from the video files with intervals of one second (total of 300 data points) were used in image analysis. Two separate training/ validation schemes were examined. In scheme A, half of all the images were used as a training set to train the model. This means that half of the images from every flare test (such as A1.1) were used for training. The other images were used to validate the model. In scheme B, for each test, all the other cases were used to train the model, which was then applied to the remaining case. This amounts to leave-one-out cross-validation, but operating on the cases rather than individual measurements. This scheme simulates the real-world application, wherein a series of different tests would be used to train the model, which would then be applied to an unknown case. The results are shown in Figures 11−14. 12650

dx.doi.org/10.1021/ie3003039 | Ind. Eng. Chem. Res. 2012, 51, 12642−12652

Industrial & Engineering Chemistry Research

Article

When training/validation scheme A is used, the model accurately predicts the combustion efficiency for both the steam-assisted and air-assisted flare. What this means is that, although the appearance of the flame varies over time, if the model is trained using some images from a given test, the other images will be similar enough (and distinct from the other tests) that the model will function properly. The different results when using scheme B (with respect to the accuracy of the model’s predictions) indicate that a wide range of conditions must be used to train the PLS model in order for it to make accurate predictions about similar (but previously unseen) flares. For flares, there were enough usable tests for the model to predict the combustion efficiency of a flare that was not in the testing set accurately.

Figure 11. Combustion efficiency (CE) predictions using the image analysis model on the air-assisted flare, with testing/validation scheme A.



CONCLUSIONS Multivariate image analysis (MIA) has been applied to flare combustion systems using both simulation and experimental testing in order to predict their resulting combustion efficiencies, as a function of crosswind velocity and steam or air injection rate. Flare images were decomposed using principal component analysis (PCA), and the resulting score vectors were used to create density score plots that describe the spectral information of the images. Segmentation techniques that use a mask in the score plot are shown to be an effective procedure way of isolating the flame region, even under complex image backgrounds. Several feature variables then were extracted from each image and then were regressed against the resulting combustion efficiency by using partial least-squares (PLS). The extracted feature variables from the flare images showed a strong correlation to the combustion efficiencies of these units. Observed and predicted values of the flare combustion efficiency at different crosswind velocities and steam or air injection rates were in good agreement. Hence, it was shown that multivariate regression models based on flare images can be used to predict the flare performance at different crosswind velocities and with different operating conditions. This work suggests that simple two-dimensional (2-D) color images of industrial flares may be a fast, accurate, and inexpensive approach for online monitoring of these complex combustion systems. This would allow for developing effective flare control and mitigation strategies. Future research work should could be made to consider the fuel flow rate, the fuel heat content, and the crosswind velocity as new parameters to improve the PLS regression.

Figure 12. Combustion efficiency (CE) predictions using the image analysis model on the air-assisted flare, with testing/validation scheme B.

Figure 13. Combustion efficiency (CE) predictions using the image analysis model on the steam-assisted flare, with testing/validation scheme A.



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Present Addresses †

Now at Quantum Reservoir Impact, Houston, TX 77010, USA. ‡ Now at Carnegie Mellon University, Pittsburgh, PA 15213, USA. ‡

Notes

The authors declare no competing financial interest.



Figure 14. Combustion efficiency (CE) predictions using the image analysis model on the steam-assisted flare, with testing/validation scheme B. (There were not enough usable tests for this approach to work with the steam-assisted flare.)

ACKNOWLEDGMENTS The preparation of this paper is based on work supported by the State of Texas through the Air Quality Research Program 12651

dx.doi.org/10.1021/ie3003039 | Ind. Eng. Chem. Res. 2012, 51, 12642−12652

Industrial & Engineering Chemistry Research

Article

administered by the University of Texas at Austin by means of a grant from the Texas Commission on Environmental Quality.



REFERENCES

(1) Shimoda, M.; Sugano, A.; Kimura, T.; Watanabe, Y.; Ishiyama, K. Prediction methods of unburnt carbon for coal fired utility boiler using image processing technique of combustion flame. IEEE Trans. Energy Convers. 1990, 5 (4), 640. (2) Yamaguchi, T.; Grattan, K. T. V.; Uchiyama, H.; Yamada, T. A practical fiber optic air-ratio sensor operating by flame color detection. Rev. Sci. Instrum. 1997, 68 (1), 197. (3) Huang, Y.; Yan, Y.; Lu, G.; Reed, A. On-line flicker measurement of gaseous flames by image processing and spectral analysis. Meas. Sci Technol. 1999, 10, 726. (4) Lu, G.; Yan, Y.; Ward, D. D. Advanced monitoring, characterization and evaluation of gas fired flames in a utility boiler. J. Inst. Energy 2000, 73, 43. (5) Wang, F.; Wang, X. J.; Ma, Z. Y.; Yan, J. H.; Chi, Y.; Wei, C. Y.; Ni, M. J.; Cen, K. F. The research on the estimation for the NOx missive concentration of the pulverized coal boiler by the flame image processing technique. Fuel 2002, 81, 2113. (6) Yu, H.; MacGregor, J. F. Monitoring flames in an industrial boiler using multivariate image analysis. AIChE J. 2004, 50 (7), 1474. (7) Szatvanyi, G.; Duchesne, C.; Bartolacci, G. Multivariate image analysis of flames for product quality and combustion control in kiln. Ind. Eng. Chem. Res. 2006, 45, 4706. (8) Kuipers, E. W.; Jarvis, B.; Bullman, S. J.; Cook, D. K.; McHugh, D. R. Combustion efficiency of natural gas flares: Effect of wind speed, flow rate and pilots. Internal Report, Shell Research and Technology Thornton and British Gas Research Centre, 1996. (9) Bandaru, R. V.; Turns, S. R. Turbulent jet flames in a crossflow: Effects of some jet, crossflow, and pilot-flame parameters on emissions. Combust. Flame 2000, 121, 137. (10) Johnson, M. R.; Kostiuk, L. W. Efficiencies of low-momentum jet diffusion flames in crosswinds. Combust. Flame 2000, 123, 189. (11) Castiñeira, D.; Edgar, T. F. CFD for simulation of crosswind on the efficiency of high momentum jet turbulent combustion flames. J. Environ. Eng. 2008, 134 (No. 7), 561. (12) Geladi, P.; Grahn, H. Multivariate Image Analysis, Wiley: Chichester, U.K., 1996. (13) Wold, S. Cross-validatory estimation of the number of components in factor and principal component models. Technometrics 1978, 4, 397. (14) Jackson, J. E. A User’s Guide to Principal Components. Wiley: New York, 1991. (15) Tobias, R. An introduction to partial least squares regression, Report TS-509; SAS Institute, Inc.: Cary, NC, April 1997. (16) Rosipal, R.; Kramer, N. Overview and recent advances in partial least squares. In Subspace, Latent Structure, and Feature Selection Techniques; Saunders, C. et al., Eds.; Lecture Notes in Computer Science, Vol. 34; Springer: New York, 2006. (17) Wold, H. Path models with latent variables: The NIPALS approach. In Quantitative Sociology: International Perspectives on Mathematical and Statistical Modeling; Blalock, H. M. et al. Ed.; Academic Press: New York, 1975. (18) Castiñeira, D.; Edgar, T. F. CFD for simulation of steam assisted and air assisted flare combustion systems. Energy Fuels 2006, 20, 1044. (19) Geladi, P.; Isaksson, H.; Lindqvist, L.; Wold, S.; Esbensen, K. Principal component analysis of multivariate images. Chemometr. Int. Lab. Syst. 1989, 5, 209. (20) Esbensen, K.; Geladi, P. Strategy of multivariate image analysis (MIA). Chemometr. Intell. Lab. Syst. 1989, 7, 67. (21) Liu, J. J.; Bharati, M. H.; Dunn, K. G.; MacGregor, J. F. Automatic masking in multivariate image analysis using support vector machines. Chemometr. Intell. Lab. Syst. 2005, 79, 42.

12652

dx.doi.org/10.1021/ie3003039 | Ind. Eng. Chem. Res. 2012, 51, 12642−12652