UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

Paper machine data analysis using wavelets Nesic, Zoran 1996

You don't seem to have a PDF reader installed, try download the pdf

Item Metadata

Download

Media
[if-you-see-this-DO-NOT-CLICK]
ubc_1996-0144.pdf [ 15.7MB ]
[if-you-see-this-DO-NOT-CLICK]
[if-you-see-this-DO-NOT-CLICK]
Metadata
JSON: 1.0065195.json
JSON-LD: 1.0065195+ld.json
RDF/XML (Pretty): 1.0065195.xml
RDF/JSON: 1.0065195+rdf.json
Turtle: 1.0065195+rdf-turtle.txt
N-Triples: 1.0065195+rdf-ntriples.txt
Original Record: 1.0065195 +original-record.json
Full Text
1.0065195.txt
Citation
1.0065195.ris

Full Text

Paper Machine Data Analysis Using Wavelets by Zoran Nesic B.Sc. University of Belgrade, Yugoslavia, 1987  A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF APPLIED SCIENCE  in THE FACULTY OF GRADUATE STUDIES DEPARTMENT OF ELECTRICAL ENGINEERING  We accept this thesis as conforming to the required standard  THE UNIVERSITY OF BRITISH COLUMBIA February 1996 © Zoran Nesic, 1996  In presenting this thesis in partial fulfillment of the requirements of an advanced degree at the University of British Columbia, I agree that the library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the Head of my Department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission.  Department of Electrical Engineering  The University of British Columbia 2356 Main Mall Vancouver, B C V6T 1Z4 Canada Date: February 1996  Abstract  The thesis describes a new approach to paper machine process data analysis using onedimensional and two-dimensional discrete wavelet transforms. These techniques have been adapted from a general theory that has been developed in recent years on the application of wavelets to signal analysis. Application areas in which the theory was first applied have included image processing and bandwidth compression for communications. Two main applications of the discrete wavelet transform have been analyzed in this thesis. First, an analysis of the use of wavelets for processing scanned data representing basis weight and moisture variations on a paper machine has been carried out. It has been shown that wavelets are effective for the detection of process signals in noisy data, so leading to better estimation and visualization of the machine direction and cross machine variations in process data. The second main application of the method has been to allow significant compression of the process data without diminishing the ability to reconstruct accurate profiles. It has been shown that the compression method can be embedded into the estimation algorithm, producing excellent results without a major expense in computation time. It has been shown that, in both applications, the new methods produce results superior to the industrially accepted procedures.  For appropriate choice of wavelets, profile estimates are  improved over those obtained using exponential filtering or other standard analysis methods. The data compression technique presents a new concept in paper machine data analysis and the author is not aware of any previous references to this subject. The ability to reduce data storage requirements is of importance in mill-wide process monitoring systems. A comprehensive analysis of the proposed algorithms has been carried out on a variety of simulated data sets for which the true process variations are known. Industrial data has also been analyzed and it is apparent that the method had many desirable characteristics.  I  Table of Contents  Abstract  II  List of Tables  Y  \  List of Figures  VI  Acknowledgments  ,  1 Introduction  '*  1  1.1  Measurement of Paper Properties  1  1.2  Review of Standard Estimation Methods  2  1.3  Motivation  3  1.4  Outline of Thesis  4  2 Wavelet Theory  5  2.1  Definitions  2.2  Wavelets as Basis Functions  6  2.3  Wavelet Properties  8  2.4  Implementation of the Discrete Wavelet Transform  9  2.5  Detection of Signals in Noise and Data Compression  11  2.5.1  12  :  De-Noising  2.5.1.1 Thresholding 2.5.1.2 Threshold  Selection  2.5.1.3 De-Noising 2.5.2  methods  of Two-Dimensional Signals  Compression  3 Analysis of Paper Machine Process Data 3.1  5  13 14 15 16  19  Matlab™ Toolbox for Wavelet Analysis of Paper Machine Data  19  3.1.1  21  Padding of the Two-Dimensional Signals  3.2  Applying Wavelet De-noising to Simulated Scanned Data  21  3.2.1  Edge Effects and Step Changes  21  Smooth Surfaces  29  3.2.3  Grade Profile Change  33  3.2.4  Data Created by the Moisture Variation Model ( M V M )  38  3.2.5  Bumps and Streaks  42  . 3.2.2  3.3  4  5  Wavelet Analysis of Industrial Data  47  3.3.1  Choice of an Estimator  47  3.3.2  Pulp Machine Moisture Content Data  47  3.3.3  Grade Change for a Paper Machine  56  Other Applications of Wavelets in Paper Making Process  59  4.1  Compression of Paper Machine Data  59  4.2  Wavelet Estimation for Control Purposes  73  Conclusion  75  References  78  A  Exponential Multiple Scan Trending (EXPO) Algorithm  82  B  Estimation and Identification of Basis Weight and Moisture Content (EIBMC) Algorithm83 B.l  B.2  Process Models  83  B . l . l Moisture  83  B.1.2  84  Basis Weight  The E I B M C Algorithm  85  B.2.1  Moisture  85  B.2.2  Basis Weight  86  References  87  C  88  Wavelet Filter Coefficients  IV  List of Tables  3.1 Edge effects: thresholds and estimated noise at each resolution level  25  3.2 Edge effects: summary of the results  26  3.3 Edge effects: best filter performances  27  3.4 Smooth surfaces: thresholds and estimated noise at each resolution level  29  3.5 Smooth surfaces: summary of the results  29  3.6 Smooth surfaces: best filter performances  30  3.7 Grade change: thresholds and estimated noise at each resolution level  34  3.8 Grade change: summary of the results  34  3.9 Grade change: best filter performances  . 34  3.10 M V M test data:, thresholds and estimated noise at each resolution level  39  3.11 M V M test data: summary of the results  39  3.12 M V M test data: best filter performances  40  3.13 Bumps and streaks: summary of the results  43  3.14 Bumps and streaks: best filter performances  44  3.15 Filters used to analyze the industrial data  48  4.16 Key table  59  4.17 Edge effects: summary of compression results  60  4.18 Smooth surfaces: summary of compression results .  60  4.19 Grade change: summary of compression results  . 61  4.20 M V M test data: summary of compression results .  61  4.21 Bump data: summary of compression results  62  4.22 Pulp machine moisture data: summary of compression results  69  7  List of Figures  1.1 Paper Machine  1  1.2 Sensor path  2  2.1 Block diagram of one-dimensional DWT.  9  2.2 Block diagram of one-dimensional inverse D W T 2.3 Block diagram of two-dimensional D W T  10 . 10  2.4 Signal in time-domain and its wavelet coefficients  11  2.5 Noisy signal in time-domain and its wavelet coefficients  12  2.6 Block diagram of a compressor/decompressor without a quantizer  16  2.7 Block diagram of a compressor/decompressor with a quantizer  17  2.8 Quantizer mapping function  17  3.1 Wavelet filters used in the analysis  20  3.2 Edge effects: CD+MD before and after the noise has been added  23  3.3 Edge effects: estimated profiles.  24  3.4 Edge effects:, estimation errors as functions of time (scan number)  25  3.5 Edge effects: estimation errors as functions of time (scan number) for a smooth estimator. .  27  3.6 Edge effects: estimation results for a smooth estimator  28  3.7 Smooth surfaces: composite profiles with and without noise. .  31  3.8 Smooth surfaces: estimated profiles  32  3.9 Smooth surfaces: estimation errors  33  3.10 Grade change data set: CD+MD before and after the noise has been added  35  3.11 Grade change test data set: estimated profiles  36  3.12 Grade change test data set: estimated profiles (images)  37  3.13 Grade change test data set: estimation errors. . .  38  3.14 M V M test data: estimation errors  40  3.15 M V M test data: estimated profiles  . 41  3.16 M V M test data: C D profiles  42  3.17 Bumps and streaks: estimated profiles  45  3.18 Bumps and streaks: estimated profiles (images)  46  3.19 Pulp machine moisture data: estimated profiles  49  3.20 Pulp machine moisture data: estimated profiles [1]  50  3.21 Pulp machine moisture data: estimated profiles [1] (images)  51  3.22 Pulp machine moisture data: estimated profiles [2]  52  3.23 Pulp machine moisture data: estimated profiles [2] (images)  53  3.24 Pulp machine moisture data: estimated profiles [3] . . .  .54  3.25 Pulp machine moisture data: estimated profiles [3] (images)  55  3.26 Paper machine basis weight data  57  3.27 Paper machine basis weight data: estimated profiles.  58  4.1 M V M data set: energy distribution  63  4.2 M V M data set: distribution of nonzero wavelet coefficients  64  4.3 2D DWT: coefficient storage  64  4.4 M V M data set: compressed data  66  4.5 M V M data set: compressed data (images)  67  4.6 Pulp machine moisture data: energy distribution  69  4.7 Pulp machine moisture data: distribution of nonzero wavelet c o e f f i c i e n t s . . . . . . . . . . . . 70 4.8 Pulp machine moisture data: compressed data \7il  71  4.9 Pulp machine moisture data: compressed data (images)  72  4.10 Smooth surfaces: high resolution M D estimation  74  4.11 Pulp machine moisture data: high resolution M D estimation  74  B . l Data indexing  83  Vlll  Acknowledgments  I thank my thesis supervisor Professor Michael S. Davies for his support and guidance. I would also like to thank Professor Guy A . Dumont for his assistance and input during the course of this work Thanks to all my friends at Pulp and Paper Centre and at the Department of Electrical Engineering for making my student days enjoyable, in particular, Anil Tuladhar, Ian McLeod, Gino Carrese, Leo Stocco, Lahoucine Ettaleb, Othman Taha, and Roger Shirt. Thanks to Brian McMillan for making the computers work. Special thanks to Jahan Ghofraniha and Scott Morgan for their friendship and, of course, for sharing their rich libraries with me. This thesis would have never been possible without the support and understanding of my friend and boss Professor Andy Black. I am grateful to my wonderful wife, Maja Krzic, for her encouragement and patience during all those long years.  \X  Chapter  1:  Introduction  Chapter 1 Introduction  1.1 Measurement of Paper Properties  In the final stage of paper manufacturing, after pulping and bleaching have taken place, the mixture of water and fibre (stock) is delivered to the paper machine where moisture is first removed by drainage and mechanical pressing and then dried to produce a sheet of paper formed at the reel. A simplified diagram of a paper machine is given in Figure 1.1.  Reel Figure 1.1: Paper Machine  To achieve uniformity and high quality in the paper produced, hundreds of control loops are needed as well as experienced operators. The most important control loops are those that maintain paper properties such as basis weight (mass offibreper unit area) and moisture content (percent of water in the overall mass of paper sheet) close to the target values assigned for the grade of paper produced. 1  Chapter  Paper Sheet  I:  Introduction  Sensor path relative to sheet Cross direction  Machine direction  Sensor off-sheet Figure 1.2: Sensor path  The sensor (see Figure 1:1) travels across the sheet (cross direction or CD). During one scan (one pass over the sheet) the sensor makes between 20 and 3000 uniformly spaced measurements (Morgan [30]). A scan period can vary from 10 to 60 seconds, depending on the installation. The speed of paper sheet is up to 25 m/s or more than a hundred times the sensor speed. The combination of the paper trajectory (in machine direction or M D ) and the sensor movement (in CD) produces the zigzag pattern of measurements presented in Figure 1.2. This non-uniformly spaced sampling represents a major problem in estimation and control. The M D variations, introduced by pressure fluctuations and consistency variations in the approach system and headbox are considered to be time-dependent, fast and independent of the C D position. Control action in M D direction is normally taken at least every scan. Changes in C D profiles, coming from nonuniformity in the headbox, slice lip, pressing and drying, are considered to be relatively slow, indeed nearly time-invariant. In practice, a C D control action is taken every 2 to 4 scans (Taylor [36], Jonsson [22]). Faster detection of M D upsets is increasingly desirable with changes in paper machine design. It is obvious (see Figure 1.2) that the measured sequence of values contains information about both M D and C D variations. Decoupling the machine direction variations from the cross direction variations is not a trivial task.  1.2 Review of Standard Estimation Methods Many estimation methods have been developed for processing scanned data since automatic 2  Chapter 1: Introduction  control of paper machines was introduced. Most industrial estimators are based on the Exponential Multiple Scan Trending (EXPO) algorithm (Dahlin [4]). Assuming slowly changing, zero mean C D profiles, this algorithm (see Appendix A ) carries out exponential filtering of M D data at each C D position and produces C D profiles that are considered reasonably acceptable by the industry. The M D estimate, defined as the mean value of a scan, is available at the end of each scan. Most recent algorithms (Lindeborg [25], Chen [2], Natarajan [32], Dumont [15], Jonsson [22], Morgan [30], Wang [37]) use stochastic models for basis weight and moisture variations. They can estimate C D (Least Squares filter) and/or M D profiles (Kalman filter), with the M D estimations updated at each data point [30, 37], leading to improved M D control bandwidth. The results of the estimation method described in this thesis are compared to the results obtained by the E X P O algorithm and the Estimation and Identification of Basis Weight and Moisture Content (EIBMC) algorithm [37, 30]. The details of the latter are given in the Appendix B .  1.3 Motivation  Recent trends in actuator design have allowed control systems for paper machines to involve an increased number of measurements and control signals. Increases in resolution for both sensors and actuators are the general trend. The increased resolution potentially leads to an improved overall control of paper machines i f the multivariable control issues can be addressed effectively [36]. The increased amount of information causes problems in data processing (slows down the estimators and controllers), visualization of the process information, and in data storage. While the advances in computer technology tend to provide us with enough computing power, the visualization and data storage demands are becoming increasingly important. The machine operators benefit from a better presentation of the process data in order to be able to detect any change in the quality of final product. A n 'enhanced' image of the two-dimensional (CDxMD) data is an important tool as part of the operator interface. The enhancements may include data filtering to provide 'smooth' images or to extract some other features, for example sudden changes in paper quality, or streaks. This type of data processing will lead to better diagnostics and feature recognition, first by the human operators and then, possibly, by the automated systems. 3  Chapter  1:  Introduction  Quality control and production monitoring is increasingly important with the introduction of new standards (ISO-9001). Being able to record production and to trace back the events that caused a change in the paper quality has become essential for the paper companies. Data compression methods that provide an efficient data storage and retrieval method are desirable. In this thesis a novel way of handling scanned paper machine data is introduced. A n effort has been made to provide a framework that unites the estimation, visualizations, and data storage of the paper machine process data using wavelet analysis techniques. These techniques have been taken from the image processing literature and adopted to our purposes.  1.4 Outline of Thesis A n introduction to wavelet theory and its applications is given in the second Chapter. The third Chapter presents the results achieved by using wavelet filtering on simulated and industrial data. Applications of wavelets in paper machine process data compression and in the control of paper machines are presented in Chapter 4. The conclusions and the final remarks are given in Chapter 5. Models of basis weight and moisture variations, as well as the E I B M C and E X P O algorithms are given in appendixes.  4  Chapter 2:  Wavelet Theory  Chapter 2 Wavelet Theory In recent years wavelet theory has been established as a unified framework for a number of techniques used in different disciplines of science: •  Signal Processing (multiresolution analysis, subband coding, pyramid schemes)  •  Applied Mathematics (multigrids for solving partial differential equations, approximation theory, grids)  •  Physics (the theory of coherent states) Statistics (regression)  The roots of what is today known as wavelet theory can be found in the work done in the early 1930's (Littlewood-Paley techniques) or even earlier (1910, work done by A . Haar [17]). In the early 1980's, these different techniques were unified as wavelet theory. At that time, the word wavelet was suggested by Alex Grossman and Jean Morlet. In mid 1980's, a new orthogonal wavelet expansion was constructed by Pierre-Gilles Lemarie and Yves Meyer [24] and, by the end of the decade, Ingrid Daubechies [5] gave a method for the construction of wavelets, nonzero only on a finite interval and with arbitrarily high but fixed regularity. Most of the work carried out in applying wavelet theory in the detection of signals in noise, and in data compression has been published in the last five years. In this Chapter an introduction to wavelet theory, one- and two-dimensional discrete wavelet transforms, wavelet filtering, and data compression is given, with an emphasis on the implementation issues.  2.1 Definitions The inner product of two elements x, y of a real or complex linear space Y is denoted (x, y) and has the following properties: (i)  (x, y) > 0 , whenever x ^ 0  (ii) (x,y) = (y,x) (iii) (ax, y) — a(x,y) , a is a scalar (iv) (x + y,z) = (x, z) + (y, z)  5  Chapter 2:  The norm, denoted ||y||, can be defined as ||y|| = (y, y) . ¥  x, y £ Y are defined to be orthogonal (denoted x _L y) if  Wavelet Theory  In an inner product space, two elements y) = 0. Let {y,-} be a set of mutually  orthogonal ((yi,yj) — 0, i ^ j) elements in an inner product space Y. Define e,; = yi/||y»:||, then each e, has unit norm and the set {e -} is called an orthonormal set. s  A Hilbert space is a linear space on which an inner product has been defined, with a norm derived from the inner product and which is complete in this norm. 1? is a Hilbert space with inner product  (x,y) = Jx{t)y(t)dt  (2.1)  where the functions x, y are defined and integrable on a domain I. I is a Hilbert space with inner 2  product  (z> v) =i=iX] oo  The condition number k of an operator T is defined as k = ||T|| • if T* =  (-)  XiVi  22  1|. A n operator is unitary  T~ . l  2.2 Wavelets as Basis Functions Any signal in a Hilbert space can be approximated by a weighted sum of basis functions. The accuracy of the approximation will depend upon the class of signals, the choice of basis functions and the number of terms. N  f(x)^^2c^i(x) Different sets of basis functions will lead to different sets of coefficients c;.  (2.3) Two classes  of representation are, for example, the sampling function and the sinusoid function. The signal representation by using sinusoids as the basis function is known as Fourier transform. This transform reveals information only about the signal's frequency content over the entire interval. On the other hand, the signal representation by a weighted sum of sampling functions will contain information only about the time domain behavior. This can be explained by comparing the functions' support (support is the interval over which a function is nonzero). The sampling function has infinitesimally small support and so produces precise time information but no frequency data. The sinusoids have infinite support and so provide spectral information but cannot establish time-domain variations. In 6  Chapter 2:  Wavelet Theory  many cases a set of basis function that is a compromise between these two extremes can be used to track both spectral and time-domain disturbances. These functions should have finite support (good time resolution) of adjustable width [18] (good frequency resolution). One such a set of basis functions is called wavelets.  Wavelets are created by scaling and  translating the same prototype function ^(x) known as the mother wavelet.  The scale factor is  normally chosen to be a power of two, yielding the desired cascade octave band-pass filter structure. For the discrete wavelet transform, all the integer shifts of *&(x) have to be considered. Therefore, the wavelet decomposition of the signal is given by (2.4) where V (x) = 2iV(2 x-k)  =2*V(2 (x-2- k))  j  j  jk  j  (2.5)  The coefficient 2~ k is a dyadic point. The multiplier 2 s is needed to make the basis orthonormal. 3  The wavelet coefficients Cjk are computed by the continuous wavelet transform, which is the inner product of the signal f(x) and the basis function &jk(x). Let the time and the frequency resolutions be defined (in rms sense) as  -OO  At =  (2.6)  oo  \  ; w)\ dt 2  J uj \G{uj)\ du 2  2  Au  (2.7)  / \G{u)\ du 2  where g(t) and G(u) are basis function and its Fourier transform respectively. The denominators represent the energies of g(t) and G(u) respectively. The product of time and frequency resolutions is lower bounded [34] (Eq. 2.8).  At Aw > ~ 2  (2.8)  Eq. 2.8 is sometimes referred to as the uncertainty principle, or Heisenberg inequality. It implies that time and frequency resolutions cannot be controlled independently.  Chapter 2:  Wavelet Theory  2.3 Wavelet Properties Up to this moment, the choice of the mother wavelet of the function ^(x)  has not been discussed. The choice  depends on the wavelet properties that are best suited to a particular problem.  Some of the most important properties of wavelets are given below. For a more detailed discussion of wavelet properties see the work of Jawerth and Sweldens [20].  Orthogonality  - /)) = <5/) directly links  the L norm of a function to the norm of its 2  wavelet coefficients by (2.9)  The fast wavelet transform, using orthogonal wavelets, is a unitary transformation. Orthogonality is also important in numerical calculations. A n error present in the initial data will not grow under the transformation, and stable numerical computations are possible. In the multiresolution analysis orthogonality means that the projection operators into the different subspaces yield optimal approximations in the L  2  Compact support is  sense. needed to guarantee that the summations in the fast W T are finite (i.e.  wavelet filters have finite impulse responses).  Symmetry  is required for the generalized linear phase filters (to prevent phase distortion).  Smoothness is  a very important property for wavelet filtering (de-noising) and compression  applications. A function is said to have M degrees of smoothness i f its M  t  n  derivative is continuous  at all points. The reconstructed version of a function will be much smoother i f the wavelet itself has a higher degree of smoothness. A higher degree of smoothness corresponds to better frequency localization of the filters.  Number of vanishing moments TY for a sequence {gk}  is given by Eq. 2.10 and is connected to  the smoothness of the wavelet (and vice versa). A higher number of vanishing moments corresponds to a higher degree of smoothness. (2.10) k  There is a trade-off among these properties and they are held to varying degrees by different families of wavelets. There are many researchers working on the wavelet design. Fundamental  8  Chapter 2:  Wavelet Theory  and pioneering work has been done by Daubechies [6, 8, 5, 7], Meyer [29, 28], Cohen [3] and Wickerhauser [38]. Having chosen the appropriate wavelet, one can move towards the implementation of the discrete wavelet transform.  2.4 Implementation of the Discrete Wavelet Transform The paper machine sensor provides a sequence of discrete measured values.  Therefore, an  efficient algorithm for the discrete wavelet transform (DWT) is needed. Mallat's pyramid algorithm [26] has been used in this thesis and it will be described in this section. Let the sequence Sk be  HP j-2 ' HP  i=  LP  i:  1=2  • ••  Figure 2.1: Block diagram of one-dimensional DWT. of the length n = 2 . N  The block diagram of a one-dimensional discrete wavelet transform (DWT)  is given in Figure 2.1. It requires two properly designed wavelet filters: G — high pass filter and H — low pass filter, both of the length M. Examples of such filters are given in Figure 3.1 and Appendix C. These two filters are not independent since:  9n  where g  n  and h  n  =  h\-n  n  (2.H)  are filter coefficients of G and H respectively. The D W T is implemented by  filtering the signal followed by down-sampling by a factor of 2 (Figure 2.1 and Eq. 2.12). M-l  HPk = ^ 9iS2k-i,  (2.12)  M-l  LPk =  ^2 i=0  hiS2k-i  If the length of the original sequence Sk is n then the lowpass (LP) and the highpass (HP) component are each of the length n/2. This procedure of filtering and down-sampling can be repeated recursively  Chapter  2:  Wavelet Theory  on the L P components of the signal. After N iterations the length of the lowpass and highpass components becomes equal to 1. The wavelet coefficients are given by:  w = {LP , HP , j=N  HP -i,  j=N  j=N  • • •, HP }  (2.13)  j=1  HP,  LPM  L P *  ,-^f2)-»>  •(t2)-^  H  H  Figure 2.2: Block diagram of one-dimensional inverse DWT. The inverse D W T is depicted in Figure 2.2. It is accomplished by upsampling of the highpass and lowpass components by inserting a zero after each element, and convolving the upsampled lowpass and highpass signals with the filters H and G respectively. The filter pair (H, (?) is related to the forward D W T filter pair (H, G) 9n=  9-n  (2-14) h  n  — h—jj  The fast D W T requires 0(n) operations. A n actual algorithm in pseudo code is given in [20]. A two-dimensional D W T can be accomplished by two separate one-dimensional transforms [18] as presented in Figure 2.3. r—  H(y )  12 along y— * -  G(y )  12 along y —  H(y  12 along y  k  H(XK)  —»>  12 along x — * -  — * •  k  — * •  k  k  —»*- i2  along x  k—  k  k  •  k  » -  G(y  k  12 along y  k  Figure 2.3: Block diagram of two-dimensional DWT. When implementing the filtering of a signal by the filters H,G a choice has to be made regarding the signal padding at the ends of the sequence. A signal padding assumption is needed because both, 10  Chapter  Wavelet Theory  2:  the filter and the signal lengths are finite. When the filter passes the last signal value (s ) n  h\ ...  S „ _ 4 S „ _ 3 S -2 n  h-2  S -l n  /13 fi4 S X\ n  £3  X2  (2.15)  ...  some values {x\, X2, £ 3 , . . . } have to be added at the end of the signal. The simplest solution would be to pad zeros. Padding zeros introduces a discontinuity at that point and it is quite unlikely that zeros are a natural extensions of the signal. Knowing how sensitive the D W T is in registering discontinuities, which is a reason for using the D W T in the first place, it is likely that this artificial insertion will introduce a significant error. x'% = S2,  A better approach is to consider the function to be periodic (i.e. x\ = s\, Again, it will cause no error only i f the behavior of the signal at s i and s  n  ...).  do match. Other much  more advanced approaches to this problem, are summarized in [20] (Meyer's boundary wavelets, Dyadic boundary wavelets).  2.5 Detection of Signals in Noise and Data Compression The signal representation in wavelet domain can become very compact with only a few coefficients being needed to represent a complex signal (Figure 2.4). B y contrast, broadband components associated with noise have their energy dispersed over the coefficients (Figure 2.5). Wavelet coefficients  Function Blocks (MakeSignal.m, WaveLab) 0  61  -1 -2 -3 -4 -5 -6 -7 -8 0.2  0.4 0.6 Time  -9  0.8  0.2  0.4 0.6 Time  Figure 2.4: Signal in time-domain and its wavelet coefficients 11  0.8  Chapter 2:  Wavelet Theory  This property makes the wavelet transform very suitable for estimation and data compression. By reducing (shrinking) the wavelet coefficients it is possible to remove the "less significant" features of a signal, most likely noise. Therefore an estimated, noiseless, signal can be obtained after the inverse D W T . The fact that the signal energy is concentrated in a relatively small percentage of its wavelet coefficients, naturally led to the implementation of wavelets in data compression. 2.5.1 De-Noising The term De-Noising describes, in an informal way, the various schemes which attempt to reject noise by damping or thresholding in the wavelet domain [11]. Different methods of wavelet filtering are very well documented in the literature. Most of the implementations are based on the work done by Donoho and Johnstone at Stanford University. They have shown [13, 12, 14, 10] that for the signal model: i — 1,..., n  Vi = fi + crei,  (2.16)  where e; ~ N(0,1) is Gaussian white noise, the non-linear wavelet-coefficient thresholding produces l  estimates / that are nearly optimal in the I sense. 2  Wavelet coefficients in presence of noise 0  Function Blocks with additive white noise 7.  -1 -2 -3 _i_  -4 -5  I  -6  1  T I'  T"  r~r  •—TA  J  A  ,,1. ,  -7  -Ir-  -8 0  0.2  0.4 0.6 Time  -9  Oi  0  0.2  0.4 0.6 Time  Figure 2.5: Noisy signal in time-domain and its wavelet coefficients 12  Oi  Chapter  2:  Wavelet Theory  The wavelet regression estimator that is suggested by the Stanford group is based on a three-step procedure: 1.  Find D W T of y =  (yi,...,y„) u = Wy  (2.17)  where operator W represents the wavelet transform; 2.  Create a new set of wavelet coefficients UJ* by a non-linear modification of UJ;  3.  Reconstruct an estimate of / b y applying the inverse D W T to UJ*: f = W u* T  (2.18)  The step 2 is the most important part of the regression. Donoho and Johnstone have proposed thresholding of wavelet coefficients as the primary method of the modification. They have established that this method works because: 1.  The D W T of Gaussian white noise is again white noise. Thus, it is evenly spread over all coefficients (every empirical wavelet coefficient contributes noise of variance a  2  but  only a very few wavelet coefficients contribute signal [14]. 2.  The real signal / is represented by a small number of nonzero coefficients.  To achieve successful estimation, a threshold has to be found based on the estimated noise level. 2.5.1.1 Thresholding methods There are two main types of thresholding: Hard Thresholding and Soft Thresholding. Hard thresholding keeps all the coefficients that are greater than the threshold and sets the others to zero. Expressed in the standard Matlab™ notation, it is given by ui*  hard  — ui • (abs(ui) > threshold)  (2-19)  Soft thresholding sets all the coefficients smaller than the threshold to zero and shrinks the others by the threshold value. It is given by ui* f = sign(u>) • (abs(uj) — threshold) • (abs(uj) > threshold) so  t  13  (2.20)  Chapter  2:  Wavelet Theory  It is very common to apply thresholding only on the higher-resolution coefficient levels, keeping the lower-frequency content of the signal intact. It is obvious that these methods are highly non-linear. The soft thresholding is the I optimal 2  non-linear function in the wavelet domain to apply i f one requires the resulting function to be at least as smooth as the original, noise-free one. Hard thresholding yields better I performance but 2  does not guarantee the smoothness property. In practice, hard thresholding may produce estimates containing "blips" and other irregularities that may not be visually appealing. Another method of thresholding, proposed by Lang et al. [23], uses a shift invariant D W T for the noise reduction. Authors have shown that this method, when combined with hard thresholding, gives better noise suppression without oversmoothing of the details. This method calculates the D W T for all the shifts of the signal, followed by the usual hard thresholding and inverse D W T . Then, it averages the resulting estimates. When properly implemented, shift invariant D W T requires O(NlogN)  operations and storage space for NlogN coefficients compared to O(N) operations and  storage space for N coefficients required for classical DWT. 2.5.1.2 Threshold Selection Donoho and Johnstone have established [13, 10] several methods of threshold selection.  Universal Threshold or VisuShrink is given by tuniv = cryj2log(n) where n is the sample size (2 ) and a is the noise-level estimate. N  (2.21) The authors have proposed  median absolute deviation (MAD) method for the noise level estimation [13] from the finest scale of empirical wavelet coefficients (ojf) a = Median(\u \)/0.67i5 f  (2.22)  This is a robust estimator which produces visually appealing (hence the name VisuShrink), very smooth, estimates.  The same threshold is applied to all the levels. It has been shown that this  method tends to underfit (oversmooth) the data [31]. SureShrink is an adaptive threshold-selection procedure [13]. A different threshold is applied at each resolution level leading to smoothness-adaptive estimates. It uses Stein's Unbiased Risk 14  Chapter  2:  Wavelet Theory  Estimate for threshold choice. A n implementation of SURE algorithm is given in the Matlab file: ValSUREThresh.m which is part of the wavelet toolbox developed at Stanford University [1]. At each level of resolution the algorithm checks whether there is enough signal to compute the SURE threshold. If there is, t  sure  is used, otherwise  t iv is used. This method also tends to underfit. un  There are some other methods, proposed by different authors. Nason [31] has suggested two methods. One uses SURE procedure on all the coefficients simultaneously. It is called GlobalSure. This procedure results in reasonable performance with respect to mean square error. The second method is based on Cross-Validation. The same paper presents a comparative study of these methods. For all the threshold-selection procedures described above, the noise model is assumed to be as in (2.16). In case of correlated noise Johnstone and Silverman [21] have proposed a level-adaptive threshold: (2.23)  tj = Sj sjl • login)  where n is the data sample size and Sj is the standard deviation of the wavelet coefficients at j  t h  level.  Another method that can be used in the case of correlated noise, i f the noise structure is known, is to use a prewhitening transformation followed by universal thresholding. However, the wavelet decomposition of the signal in the original domain may possess sparsity properties that are lost in the prewhitening transformation, and the advantages of using wavelet de-noising on the prewhitened data would then be diminished [21].  2.5.1.3 De-Noising of Two-Dimensional Signals  The same methods described in the preceding sections can be used for two-dimensional signals. Donoho [10] proposed the same three step procedure described in 2.5.1, but applied to twodimensional image data  i,...  ,mi  1,... ,m  2  where zi  iidN(0,1)  is Gaussian white noise and the sample size n is given byn = m\m2. 15  (2-24)  Chapter  2:  Wavelet Theory  2.5.2 Compression Data compression is one of the most common applications of wavelet theory. The ability of the D W T to compress the signal energy in a small number of coefficients makes it a natural choice for data compression. There are two basic types of data compression, lossless and lossy. Lossless compression is used when an exact reconstruction of the signal is needed.  In con-  trast, lossy compression allows small errors in the reconstruction, in an attempt to achieve higher compression ratios.  DWT  H Threshold H  Encoder Sparse Matrix Storage  Compressor  IDWT  Decoder  Decompressor  /  Figure 2.6: Block diagram of a compressor/decompressor without a quantizer.  In this thesis, only lossy.compression will be considered. The reason for this should be apparent from the previous section (2.5.1). The goal is to try to reconstruct the" signal from the noisy data. Most of the time it is only the reconstructed signal that is of interest and the "loss" should be the noise. The de-noising procedures reduce the number of nonzero elements, therefore making it an , excellent candidate for compression. There are two possible methods to implement data compression on the process data. One method is to use the procedure depicted in Figure 2.6. The encoder's role is to losslessly compress the sparse matrix of quantized coefficients. For the fast execution, simple run-length (Pratt [33]) coding of zero-valued coefficients has been proven very effective. This simple procedure, when applied to the industrial data used later in this thesis, reduces the percentage of nonzero coefficients in the DWT to 5% (i.e. the signal is represented by 5% of its D W T coefficients). The remaining coefficients can be stored in a Matlab™ sparse matrix achieving the compression ratio of 13:1. Because the wavelet transform and de-noising procedure are an integral part of the estimator, no additional computing time is required. 16  Chapter  2:  Wavelet Theory  Another method used in image compression is shown in Figure 2.7. CO'  DWT  H Threshold H  Quantizer  H  Encoder Sparse Matrix Storage  Compressor co'.  s(\,yj'  IDWT M Dequantizer—  Decoder  Decompressor  Figure 2.7: Block diagram of a compressor/decompressor with a quantizer.  This method adds one additional step,  quantization, into the compression procedure.  By  quantization of nonzero D W T coefficients much higher compression ratios are achievable at the expense of an additional error. A  quantizer is a function that maps many input values into a smaller  set of output values. This mapping is generally a staircase function (Figure 2.8). It takes a continuous variable  u and maps it into a discrete variable u* which takes values from a finite set of reconstruction  levels { r i , T2,  TL} (Jain [19]).  t  k  U  Figure 2.8: Quantizer mapping function.  If the value of  u is between two decision levels (tk,tk+i] the quantizer will map it into a  single reconstruction level tk (Figure 2.8). The quantizer mapping is irreversible — for the given reconstruction level (a quantizer output) the input value cannot be uniquely determined. Thus, an additional distortion of the original signal is introduced. Therefore, the quantizer should be designed in such a way that it minimizes the distortion. In the case of wavelet transform, where the significant information about the signal is carried by the high magnitude wavelet coefficients, it is obvious that the quantizer should have more reconstruction levels dedicated for those coefficients than for the low 17  Chapter  2:  Wavelet Theory  magnitude ones. One design method, the optimum mean square quantizer (LLOYD-MAX) design, is given in [19]. The decompressor (Figure 2.7) applies decoding (unique) and dequantization (non-unique) to come up with a set of wavelet coefficients CJ* (note: CJ* ^ CJ*). The compression of process data differs from the compression in image processing. The image processing compression should be optimized in a sense of accounting for the Human Visual System (Marr [27]) while the goal of process data compression is to preserve the information about the process. For additional information on wavelet data compression one can refer to the overview given by Hilton et al. [18].  18  Chapter  3:  Analysis of Paper Machine Process Data  Chapter 3 Analysis of Paper Machine Process Data The first part of this chapter will be used to introduce the application of wavelet estimation methods to paper machine data analysis. A Matlab™ toolbox has been developed so that different algorithms, described in the previous chapter, can be tested. Test data sets have been used to verify the algorithms and to compare their results to those obtained by using the industrially accepted estimators (EXPO, ED3MC).  3.1 Matlab™ Toolbox for Wavelet Analysis of Paper Machine Data A set of Matlab™ functions and scripts has been developed to provide quick access to variety of different wavelet de-noising techniques and different test data sets. Even though this software is based on the Rice University Wavelet Toolbox [16] and WaveLab [1] only three functions from those toolboxes are used in the final version of the software (two-dimensional D W T and fDWT from the Rice Wavelet Toolbox and SURE threshold selection from the WaveLab). Most of the functions have been rewritten to adjust the algorithms to handle paper machine data (two-dimensional data with x and y dimensions of different, non-power of two, lengths). B y selecting different options in the front-end interface procedures many different setups can be tried . Two different thresholding methods can be used: 1.  hard thresholding,  2.  soft thresholding.  Four different methods of noise estimation and threshold selection are available: 1.  standard deviation of the wavelet coefficients at the finest resolution multiplied by y/2 * log  (n),  2.  VisuShrink procedure (see Chapter 2)  3.  M u l t i M A D , resolution-dependent thresholding using M A D noise estimator (Eq. 2.22) to estimate the noise standard deviation at each resolution level,  4.  SureShrink (see Chapter 2)  19  Chapter  3:  Analysis of Paper Machine Process Data  and filters created from four different mother wavelets can be used: 1.  Haar wavelet,  2.  Daubeshies orthogonal wavelets,  3.  Coiflets,  4.  Symmlets.  The above wavelet filters are depicted in Figure 3.1. For more information about the mother wavelets used the reader may refer to the original papers [8, 5, 9]. Haar  Daubechies 4  Figure 3.1: Wavelet filters used in the analysis. The results of wavelet de-noising procedures are compared to the estimations obtained by the E X P O algorithm and, in some cases, to those obtained by the E I B M C algorithm. Only a few tests 20  Chapter  Analysis of Paper Machine Process Data  3:  were done using the E I B M C since, for a fair comparison, the E f B M C algorithm needs to be properly tuned and the test data sets have to be created by using the proposed models (see Appendix B). The verification process called for some extreme-case input signals. Hence, not all data sets used for the wavelet estimator testing were created by the given models. Note that both E I B M C and E X P O are recursive algorithms, whereas wavelet filtering is a batch operation at present.  3.1.1 Padding of the Two-Dimensional Signals As mentioned in the previous chapter the fast DWT and IDWT procedures need input data sets of dimensions  D(2  L l  ,  2 ) . Paper machines rarely give this kind of output. Therefore, some kind of 1,2  extrapolation is needed so data can be padded up to the proper size. The padding by the repetition of the last data point in each particular direction has been chosen because it is simple to implement and it reduces the boundary errors due to the periodic nature of the input signals assumed by the D W T algorithm. The extra data points are removed after the de-noising. This procedure caused many of the boundary wavelet coefficients to be set to zero (signal becomes constant in that direction setting high-pass signal component to zero). Thus, the original M A D algorithm (Eq. 2.22) had to be modified so that all the zero coefficients are removed before the median is calculated.  3.2 Applying Wavelet De-noising to Simulated Scanned Data Verification of the proposed algorithms has been conducted on the various sets of the simulated data. Data sets were designed in such a way that they provided signals with enough complexity so that the estimator could be properly evaluated. The test data also contained the features that commonly appear on a paper sheet (streaks, bumps and grade changes). A l l the data sets have been zigzag sampled and presented to the estimator as i f they had been coming from a paper machine scanner. The data processing has been done in batch mode. A matrix of data points of the size D(n,m) is taken, padded to the size D\ (2 , 2 ), N  M  where 2  N  >n and 2  M  > ra, and sent to the estimator.  3.2.1 Edge Effects and Step Changes The first data set has been created in such a way that the edge effects (errors from the signal periodicity assumed by DWT) and the effects of sharp (sudden) changes in profiles can be studied. The target values (noiseless data) are shown in Figure 3.3. This data set' contains a very well 21  Chapter  3:  Analysis of Paper Machine Process Data  defined C D profile which differs greatly at the edge sensors (sensors #1 and #120) and it has a sharp discontinuity at the sensor #60. The M D profile is a sine function of time, with an amplitude of 1. At the M D position 29 the C D profile flattens up by 50% (sudden change in C D profile — the amplitude being changed from 1 to 0.5). The number of scans for all test data sets is 58 and the number of sensors is 120. Consequently, four data points were padded to each end of every C D profile and 3 scans were padded before and after the real signal. This signal is corrupted with white noise Zi i  lt 2  ~ N(0,a ) 2  where a = (0.2) . The composite 2  2  signal (CD + M D ) , before and after adding the noise, is shown in Figure 3.2. This noisy signal has been processed by using different wavelet estimators and the results compared to those obtained by the E X P O . For a signal containing such abrupt changes it was expected that the Haar wavelet would perform very well. This has been confirmed by the experiments (see Tables 3.2 and 3.3). The estimated profiles as well as the real signals are given in Figure 3.3. The M D profiles are calculated as scan averages.  22  Chapter  3:  Analysis of Paper Machine Process Data  MD+CD Target Profile  0  Sensors  o  Raw MD+CD Profile  Figure 3.2: Edge effects: CD+MD before and after the noise has been added.  23  Chapter  3:  Analysis of Paper Machine Process Data  Chapter  3:  Analysis of Paper Machine Process Data  Comparison of mean scan errors  Figure 3.4: Edge effects: estimation errors as functions of time (scan number). The Haar wavelet given in Table 3.3 has been used. The thresholds and the estimated noise standard deviations, used for de-noising, are given in Table 3.1. Levels  Threshold  <r(<7 = 0.2)  1  0.75  0.19  2  0.64  0.18  3  0.63  0.20  4  0.50  0.19  5  0.25  0.12  Figure 3.4 depicts the mean  Table 3.1 Edge effects: thresholds and estimated noise at each resolution level. estimation error (per scan) as a function of time (scan position). The error calculation formula is given in Equation 3.1  CDerrij = \CDij — CDestimatedij\, i = 1... n; j = 1... m where CDerrij is the absolute value of estimation error at the i  t h  scan and j  t h  (3.1)  sensor. For each  scan the mean value has been plotted.  I CDerrj = — > CDerri, m  j=l  25.  (3.2)  Chapter  J  = ^ 2 2  C  D  2  (  I  R  R  3? Analysis of Paper Machine Process Data  ' - J  (3-3)  The overall errors (Equation 3.3) were: Jwav  — 0.86j (3-4)  Jexpo = 9.6  and the mean error at each scan was less then 1% of the signal amplitude. Therefore, the signal was successfully recovered from the noisy data. A number of tests using different methods of thresholding, threshold selections, and with a variety of different mother wavelets of different lengths, has been run for each data set. Some of the results are summarized in the Table 3.2. For the Symmlets and Coiflets the.number in the fourth column represents the wavelet type not the actual length. Two tables, summarizing the test results, are given for each data set. The first table always shows the results obtained by the same set of the wavelet filters and is used for comparison of the filter performances on different data sets. The second table displays the representatives of four different wavelet filters that have produced the best results, i.e. smallest errors, for that data set. The Haar wavelet used in this example has produced the best results, as expected, but a number of different configurations have performed very well. When very smooth wavelet filters have been used, the estimation error would tend to go up with an increased number of resolution levels. This happens when thresholding starts to affect the low-frequency part of the signal (trying to smooth out the signal made of sharp edges). Thresholding  Threshold  Wavelet  Resolution  method  estimation  length/type  levels  Daubechies  soft  std  6  • 3-  0.85  Daubechies  soft  MultiMAD  6  3  1.14  Symmlet  soft  std  6  3  0.85  Symmlet  soft  MultiMAD  6  •3  1.19  Coiflet  soft  std  3  3  0.83  Coiflet  soft  MultiMAD  3  3  1.27  Haar  soft  MultiMAD  2  2  0.46  Wavelet  Table 3.2 Edge effects: summary of the results 26  .  J J expo  '  Chapter  3:  Analysis of Paper Machine Process Data  It can be seen that the errors at the edges of the signal are very small. This may be attributed to the signal padding and to the wavelet length. For a filter length of two, only one data point has to be padded to the signal at any time — consequently only one point can be wrong. J  Thresholding  Threshold  Wavelet  Resolution  method  estimation  length/type  levels  Daubechies  hard  MAD  4  2  0.70  Symmlet  hard  MAD  6  3  0.72  Coiflet  hard  std  2  3  0.66  Haar  soft  MultiMAD  2  5  0.09  Wavelet  t  « expo.  Table 3.3 Edge effects: best filter performances. To contrast the results shown in Figures 3.3 and 3.4, in Figures 3.6 and 3.5 are shown the results of the estimation obtained by using the Coiflet filter described in Table 3.3. As expected, the filter with a higher degree, of smoothness yielded much smoother estimates causing a decrease of 'sharpness' at the points of discontinuity (sensor #60 and scan #29). The overall estimation error was higher than in the first case but still significantly better than for the E X P O estimator. The filter length has also contributed to the edge effect errors. Comparison of mean scan errors  Figure 3.5: Edge effects: estimation errors as functions of time (scan number) for a smooth estimator.  27  Chapter  28  3:  Analysis of Paper Machine Process Data  Chapter  3:  Analysis of Paper Machine Process Data  3.2.2 Smooth Surfaces The second data set, depicted in Figure 3.7, has been created with the intention of testing the wavelet estimator on some smoother surfaces with a constant signal to noise ratio. The profiles CD = 1-sin 6 D,  6 D£[0M]  C  C  (3.5)  MD = 1 • sin 0 D, M  were corrupted by the additive white noise Zi i  0 D € [0, 2TT] M  *~ 7Y(0, <J ) where 2  lt 2  a = ( 0 . 5 ) which was much 2  2  higher than in the first data set. The estimation has been done by the Symmlet filter using M u l t i M A D procedure (see Table 3.5). The thresholds and the estimated noise are given in Table 3.4. The estimation errors are shown in Figure 3.9 while the estimated profiles are given in Figure 3.8.  a(a  Levels  Threshold  1  1.56  0.40  2  1.44  0.40  3  1.87  0.60  = 0.5)  Table 3.4 Smooth surfaces: thresholds and estimated noise at each resolution level.  Thresholding  Threshold  Wavelet  Resolution  method  estimation  length/type  levels  Daubechies  soft  std  6  3  0.51  Daubechies  soft  MultiMAD  6  3  0.51  Symmlet  soft  6  3  0.39  Symmlet  soft  MultiMAD  6  3  0.39  Coiflet  soft  std  3  3  0.42  Coiflet  soft  MultiMAD  3  3  0.42  Haar  soft  MultiMAD  2  2  0.95  Wavelet  . std  J J expo  Table 3.5 Smooth surfaces: summary of the results  Even though the signal to noise ratio was low, the estimated profiles were very close to the original ones and, due to the filter smoothness, very smooth. The overall errors were: Jwav  6.22,  =  (3.6) J expo — 16.03  29  Chapter  3:  Analysis of Paper Machine Process Data  Jwav — 0.3QJ~expo  Wavelet  Thresholding  Threshold  (3-7)  Wavelet  Resolution  J J expo  method  estimation  length/type  levels  Daubechies  soft  std  6  3  0.51  Symmlet  soft  MultiMAD  6  3  0.39  Coiflet  soft  MultiMAD  3  3  0.42  Haar  soft  MultiMAD  2  2  0.95  Table 3.6 Smooth surfaces: best filter performances.  30  Chapter  3:  Analysis of Paper Machine Process Data  Figure 3.7: Smooth surfaces: composite profiles with and without noise.  31  Chapter  3:  Analysis of Paper Machine Process Data  Comparison of mean scan errors  Figure 3.9: Smooth surfaces: estimation errors.  3.2.3 Grade Profile Change  This data set, shown in Figure 3.10, has been used with the wavelet estimator to see how the estimator behaves when a step change of amplitude 2 is introduced in the M D profile. The profiles  CD  = 0.3  • cos 6 , CD  OCD € [0, 2TT] (3.8)  MD  = 0.3 • s i n 6 D , M  0MD€[O,2TT]  were corrupted by additive white noise with variance that was very high compared to the profile amplitudes: a  2  =  (0.2) . 2  The estimation has been done using the Symmlet filter and M u l t i M A D  procedure (see Table 3.9). The thresholds and the estimated noise are given in Table 3.7. It is interesting to note that the VisuShrink threshold that would apply to this case was 0.66. Table 3.7 shows that a threshold of this size was applied only at one resolution level. A l l the other threshold levels were lower. Therefore, a VisuShrink estimator would oversmooth the signal in this case. This agrees with the theory given in Chapter 2. The estimation errors are shown in Figure 3.13 while 33  Chapter  Table 3.7  3:  Analysis of Paper Machine Process Data  Levels  Threshold  <j(cr = 0.2)  1  0.61  0.16  2  0.60  0.17  3  0.66  0.21  4  0.49  0.19  5  0.38  0.19  Grade change: thresholds and estimated noise at each resolution level.  the estimated profiles are given in Figures 3.11 and 3.12. The overall error for the E X P O filter was  Jexpo — 6.42. Wavelet  Thresholding  Threshold  Wavelet  Resolution  J Jexpo  method  estimation  length/type  levels  Daubechies  soft  std  6  3  0.49  Daubechies  soft  MultiMAD  6  3  0.39  Symmlet  soft  std  6  3  0.50.  Symmlet  soft  MultiMAD  6  3  0.38  Coiflet  soft  std  3  3  0.51  Coiflet  soft  MultiMAD  3  3  0.38  Haar  soft  MultiMAD  2  2  0.71  Table 3.8 Grade change: summary of the results  Wavelet  Thresholding  Threshold  Wavelet  Resolution  j J expo  method  estimation  length/type  levels  Daubechies  soft  MultiMAD  6  4  0.20  Symmlet  soft  MultiMAD  8  5  0.17  Coiflet  soft  MultiMAD  2  4  0.19  Haar  soft  MultiMAD  2  3  0.50  Table 3.9 Grade change: best filter performances.  34  Chapter 3: Analysis of Paper Machine Process Data  MD+CD Target Profile  0  Sensors  o  Raw MD+CD Profile  0  Sensors  o  Figure 3.10: Grade change data set: CD+MD before and after the noise has been added.  35  Chapter 3: Analysis of Paper Machine Process Data  36  Chapter  S  O  O  O  O  O  O  CN  >r>  3:  Analysis of Paper Machine Process Data  —  SUBOC  ui  o  j" j]  o C A  ed -a  SUBDS  SUBDS  37  Chapter 3: Analysis of Paper Machine Process Data  Figure 3.13: Grade change test data set: estimation errors.  3.2.4 Data Created by the Moisture Variation Model (MVM) This data set (the profiles are shown in Figure 3.15) has been created by the MVM. The input parameters and the profiles are given by: CD = 1 • sin 6 D exp (-0.03 * k) C  k = 1... 58, 6D G C  (scan#)  [0, 6TT] I  u =  (3-9)  0.5,  w ~ iV(0,0.1 ) 2  The estimation has been done by the Coifletfilterusing MultiMAD procedure (see Table 3.12). The thresholds and the estimated noise are given in Table 3.10. The VisuShrink threshold that would apply to this case was 0.27. Again, only one threshold (Table 3.10) used in this case was of that size. The threshold used on the lowest resolution level was much higher than the other two and was 38  Chapter 3: Analysis of Paper Machine Process Data  the result of the estimator trying to filter out the coloured noise introduced in (3.9). The estimation Levels  Threshold  a  1  0.20  0.05  2  0.31  0.09  3  0.65  0.21  Table 3.10 M V M test data: thresholds and estimated noise at each resolution level.  errors are shown in Figure 3.14 while the estimated profiles are given in Figure 3.15. The overall errors for the E X P O and the E J M B C estimators were Jexpo  =  11.43  (3.10) Jeibmc — 0.94  • Jexpo  It can be seen, from the given results, that the wavelet estimator performed very well compared to the both E X P O and E f B M C algorithms regardless of the fact that a priori knowledge of the moisture model and coloured noise has not been used. A possible improvement would be to include the noise model in the threshold-estimation procedures by increasing or lowering the threshold, at different resolution levels, by prescheduled coefficients based on the noise model. The wavelet estimated profiles were much smoother and without the artificial C D shapes that can be seen on the E f B M C estimated results (see Figure 3.16). Therefore, i f the estimated signals were to be used for control, the smooth profiles would prevent unnecessary C D control actions.  Resolution  j  Thresholding  Threshold  Wavelet  method  estimation  length/type  levels  Daubechies  soft  std  6  3  Daubechies  soft  MultiMAD  6  3  Symmlet  soft  std  6  3  Symmlet  soft  MultiMAD  6  3  0.74  Coiflet  soft  std  3  3  0.76  Coiflet  soft  MultiMAD  3  3  0.67  Haar  soft  MultiMAD  2  2  0.89  Wavelet  Jexpo  Table 3.11 M V M test data: summary of the results 39  0.78 0.69 ;  0.79  Chapter 3: Analysis of Paper Machine Process Data  j  Thresholding  Threshold  Wavelet  Resolution  method  estimation  length/type  levels  Daubechies  soft  std  6  5  0.66  Symmlet  soft  std  6  5  0.68  Coiflet  soft  MultiMAD  3  3  0.67  Haar  soft  MAD  2  5  0.76  Wavelet  Table 3.12  M V M test data: best filter performances.  Comparison o f mean scan errors i  1  1  —  I  Exp Wav Eimc  1 I  1 I  \  I  y' s. . .  V • ^~*y  ;../.>  _  j_  11  0  :  i  10  20  Figure 3.14:  i 1 30  i  Scan  1 40  M V M test data: estimation errors.  40  i 1 50  •  ' 60  Chapter  3:  Analysis of Paper Machine Process Data  Chapter  3:  Analysis of Paper Machine Process Data  CD profile #30  Scan number Figure 3.16: M V M test data: CD profiles.  3.2.5 Bumps and Streaks The ability of the wavelet estimator to localize some characteristic patterns on a paper sheet has been tested using the data set given by Equation 3.11. This data set simulates a narrow streak in M D and two actuator bump changes. The simulation results are given in Figures 3.17 and 3.18, and performances of different estimators and their respective errors of estimation are shown in Tables 3.13 and 3.14. Results shown in the Figures are obtained by the Coiflet estimator whose parameters are given in Table 3.14.  Hard  thresholding has been used with the threshold set at 0.754 and the estimated noise standard deviation was 0.18 (true value was 0.2). This signal was a very difficult test for the estimators because of its smoothness combined with some very short and/or narrow features. If the signal had contained flat profiles with the same bumps and the streak, then the Haar wavelet would have performed the best. In this case, as in the most of the other tests, the performance of the selected group of wavelet estimators is similar to the best performances for that data set. The wavelet estimator has successfully localized the test patterns of short duration. Estimation of the streak is not as effective as by the E X P O algorithm, although the contrast between smooth surfaces, the bumps, and the streak is impressive. E X P O is effective in detecting these variations since it carries on no C D smoothing. 42  Chapter 3: Analysis of Paper Machine Process Data  0.2 • cos 0CD,  6CD e [0, 2TT]  MD = 0.3 • sin6 D,  OMD G [0,2TT]  CD =  M  CT  1  st  = 0.2  bump : amplitude = 1 scan  # 3 0 . . . 39  sensor # 3 0 . . . 34 2  nd  bump:  (3.11)  amplitude = — 1 scan  # 2 0 . . . 23  sensor # 8 0 . . . 85  Streak : amplitude = 0.6 scan # 1 . . . 58 sensor # 6 0 . . . 61  Thresholding  Threshold  Wavelet  Resolution  method  estimation  length/type  levels  Daubechies  soft  std  6  3  0.82  Daubechies  soft  MultiMAD  6  3  0.99  Symmlet  soft  std  6  3  0.86  Symmlet  soft  MultiMAD  6  3  0.99  Coiflet  soft  std  3  3  0.85  Coiflet  soft  MultiMAD  3  3  1.00  Haar  soft  MultiMAD  2  2 .  1.04  Wavelet  Table 3.13 Bumps and streaks: summary of the results  43  J Jexpo  Chapter  3:  Analysis of Paper Machine Process Data  Thresholding  Threshold  Wavelet  Resolution  method  estimation  length/type  levels  Daubechies  hard  MAD  6  3  0.74  Symmlet  hard  MAD  8  3.  0.78  Coiflet  hard  std  .2  3  0.73  Haar  hard  std  2  2  0.85  Wavelet  Table 3.14  Bumps and streaks: best filter performances.  44  J Jexpo  Chapter 3:  Analysis of Paper Machine Process Data  s  -a c c3  45  Chapter  46  3:  Analysis of Paper Machine Process Data  Chapter  3:  Analysis of Paper Machine Process Data  3.3 Wavelet Analysis of Industrial Data 3.3.1 Choice of an Estimator For almost all the tested cases, the selected group of wavelet estimators performed significantly better than the E X P O filter. The purpose of this testing was to show that even though we cannot choose a single optimal wavelet estimator for all the signals, one can be chosen from the selected group (given in the tables showing the summary results) and still achieve a good performance. When the features that have to be extracted are known (when looking for certain patterns on the sheet) an appropriate wavelet estimator can be selected from the table showing the best filters for each particular job. In this way a substantially better performance (compared to the conventional estimators) can be achieved, as has been shown for the given examples. This approach to the filter selection is now used when the estimators to be implemented on the industrial data are selected.  3.3.2 Pulp Machine Moisture Content Data This data set is collected from a pulp machine manufacturing thick pulp [22] with the target value of 640.6 g/m . Slice lip opening at the headbox was 15 cm. The headbox was 4.43 m wide with 2  design speed of 184 m/min. The C D moisture control system was a Devronizer with 25 actuators each of a width of 15.24 cm. One scan took 34 seconds to complete. It has been recorded (from the actuator positions) that two C D bumps have occurred: 1.  sensor #20, scans #19—31,  2.  sensor #87, scans #13—30.  Three different wavelet estimators were used and their results were compared to the E X P O and E I B M C estimators. The chosen estimators are given in Table 3.15. The first estimator produced the results given in Figures 3.19, 3.20 and 3.21. This filter yielded very smooth results. The M D streak that can be seen in the scan #8 was picked up by the C D profile. It was concluded that this filter oversmooths data in the low-frequency range and consequently the number of resolution levels for the other two cases have been dropped to two. It has been observed, during the experiments, that the adaptive wavelet filtering methods work more reliably when the number of resolution levels has been kept low. It is well known [13, 10] that for adaptive filtering to work well, the number  47  Chapter  3:  Analysis of Paper Machine Process Data  of resolution levels used for de-noising must be much smaller than the maximum possible number of decomposition levels. The other two filters (shown here for the comparison) have been more effective in the localization of bumps and the streaks on the paper sheet. The results obtained by the second filter have been shown in Figures 3.22 and 3.23 while Figures 3.24 and 3.25 demonstrate the results achieved by the Symmlet filter. Thresholding  Threshold  Wavelet  Resolution  method  estimation  length/type  levels  Daubechies  soft  std  6  3  Daubechies  soft  MultiMAD  6  2  Symmlet  soft  MultiMAD  6  2  Wavelet  Table 3.15 Filters used to analyze the industrial data.  The wavelet filters have shown a superior performance compared to the conventional estimators when it came to the bump localization (see Figures 3.23 and 3.25). The size of the real bumps, based on the known positions of the actuators and their widths, was 13x4.8 and 18x4.8 for the first and the second bump respectively. Both, E X P O and E f B M C have shown signs of dispersing the bumps over a much greater number of scans although the width of the bumps was preserved. The wavelet estimators have detected the bump amplitude much better than the other two estimators.  48  Chapter  3:  Analysis of Paper Machine Process Data  Raw MD+CD Profile  0  o  Sensors  MD+CD Profile estimate (Wav)  Figure 3.19: Pulp machine moisture data: estimated profiles.  49  Chapter  3:  Analysis of Paper Machine Process Data  (%) ainjsio^  (%) SJTUSIOW  50  Chapter  3:  Analysis of Paper Machine Process Data  Chapter  3:  Analysis of Paper Machine Process Data  Chapter  3:  Analysis of Paper Machine Process Data  Chapter  (%) airnsiow  3:  Analysis of Paper Machine Process Data  (%) airusiopv  54  Chapter  3:  Analysis of Paper Machine Process Data  >  i  SUBDS  SUB3S  55  Chapter  3:  Analysis of Paper Machine Process Data  3.3.3 Grade Change for a Paper Machine This data set is collected from a newsprint machine during a grade change with the dry weight control system in auto mode. Data used for the analysis was the basis weight in grams per square meter. There were 141 sensors and data was collected during 189 scans. This data set was much larger than those previously analyzed so that it was possible to select estimators that used more resolution levels for de-noising. A Coiflet3 filter has been used and soft thresholding has been applied on 5 resolution levels using a threshold t = (j-^/2 * log (n * m) where the standard deviation has been estimated from the wavelet coefficients at the finest resolution level: a = std(u>f). The results are shown in Figures 3.26 and 3.27. The estimated C D profiles look very smooth, and it is possible to see that some of the M D change has found its way into the C D profile. This has happened because the M D profiles were calculated (and updated) once per scan. A n alternative method of M D profile estimation will be proposed in the next Chapter in order to avoid this difficulty.  56  Chapter  3:  Analysis of Paper Machine Process Data  Raw MD+CD Profile  MD+CD Profile estimate (Wav)  Figure 3.26: Paper machine basis weight data.  57  Chapter  3:  Analysis of Paper Machine Process Data  Chapter 4:  Other Applications of Wavelets in Paper Making Process  Chapter 4 Other Applications of Wavelets in Paper Making Process  4.1 Compression of Paper Machine Data The theoretical background needed for the understanding of the compression algorithms used in this chapter has been given in Section 2.5.2. In this chapter some of the more practical aspects of the process-data compression method, as well as the results achieved will be given. Both compression methods given in Section 2.5.2, (compression with and without quantizer) have been used. It has been assumed that the input data has been measured and stored in 16-bit format. A l l the calculations have been done using double precision (Matlab™'s default precision).  The software has been  written in such a way that quantizer resolution could assume any value. Only three resolutions have been used for the test purposes: 8—, 12— and 16-bit. The 16-bit resolution has been used as the reference, yielding the compression ratio of 1:1. The 12—bit resolution corresponded to a compression ratio of 1.33:1 and 8—bit resolution resulted in a compression factor of 2:1. A l l the above mentioned compression ratios are in addition to the compression achieved by the thresholding of wavelet coefficients and the encoding of data into sparse matrices. Case  Wavelet  Threshold  Wavelet  Resolution  estimation  type/length  levels  1.  Daubechies  std  6  3  2.  Daubechies  MultiMAD  6  3  3.  Symmlet  std  6  3  4.  Symmlet  MultiMAD  6  3  5.  Coiflet  std  3  3  6.  Coiflet  MultiMAD  3  3  7.  Haar  MultiMAD  2  2  Table 4.16 Key table The examples used in Chapter 3 have been used to test the ability of the wavelet compressor to compress different types of signals. The results are summarized in Tables 4.17 to 4.21. The wavelet filters used here are given in Table 4.16. 59  Chapter  compress,  compress,  Jcompress Jexpo  ratio  8-bit  12-bit  16-bit Case  Other Applications of Wavelets in Paper Making Process  4:  Jcompress Jexpo  ratio  compress,  Jcompress Jexpo  ratio  1.  30.00  0.85  40.00  0.85  60.00  1.17  2.  50.80  1.15  67.74  1.15  101.61  1.43  3-  28.06  0.85  37.42  0.85  56.13  1.02  4.  52.33  1.19  69.77  1.19  104.66  1.32  5.  27.08  0.82  36.11  0.83  54.16  0.85  6.  52.73  1.27  70.30  1.27  105.45  1.29  7.  13.57  0.46  18.09  0.46  27.13  0.58  Table 4.17 Edge effects: summary of compression results  Case  compress,  compress,  Jcompress Jexpo  ratio  8-bit  12-bit  16-bit  ratio  Jcompress Jexpo  compress, ratio  Jcompress Jexpo  1.  50.43 .  0.51  67.25  0.51  100.87  0.66  2.  50.80  0.51  67.74  0.51  101.87  0.66  3.  49.36  0.39  65.82  0.39  98.72  0.51  4.  51.56  0.39  68.74  0.39  103.11  0.51  5.  49.01  0.42  65.35  0.42  98.03  0.43  6.  53.54  0.42  71.38  0.42  107.08  0.42  7.  13.57  0.95  18.09  0.95  27.13  0.96  Table 4.18 Smooth surfaces: summary of compression results  The actual quantization used here is given by this four-step procedure (for a data set D): 1. find t i m  n  = min(D) and t  max  2.  A t = (tmax ^min)/(2  3.  Dj = D —  4.  D = c  —  max(D)  l)  t in m  floor(Di/At)  where the variable bits represents the resolution in the number of bits and the function floor rounds the given number to the first integer smaller or equal to the number. 60  Chapter 4:  Other Applications of Wavelets in Paper Making Process  The reconstruction (de-quantization) algorithm is given by D' = D * At + t in c  compress, ratio  8-bit  12-bit  16-bit Case  (4-1)  m  Jcompress Jexpo  compress, ratio  Jcompress Jexpo  compress,  Jcompress Jexpo  ratio  1.  19.39  0.49  25.85  0.49  38.77  0.51  2.  50.08  0.39  67.74  0.39  101.61  0.42  3.  16.93  0.50  22.58  0.50  33.87  0.53  4.  51,56  0.38  68.74  0.38  103.11  0.42  5.  17.66  0.51  23.55  0.51  35.33  0.53  6.  53.54  0.38  71.38  0.38  107.08  0.90  7.  13.57  0.71  18.09  0.71  27.13  0.71  Table 4.19 Grade change: summary of compression results  Case  compress, ratio  8-bit  12-bit  16-bit Jcompress Jexpo  compress, ratio  Jcompress Jexpo  compress, ratio  J compress Jexpo  1.  6.98  0.78  9.31  0.78  13.97  0.82  2.  48.00  0.69  64.00  0.69  96.00  0.74  3.  6.87  0.79  9.16  0..79  13.74  0.85  4.  44.38  0.74  59.18  0.74  88.77  0.80  5.  6.81  0.76  9.09  0.76  13.63  0.76  6.  45.96  0.67  61.28  0.67  91.91  0.67  7.  12.58  0.89  16.78  0.89  25.17  0.89  Table 4.20 M V M test data: summary of compression results . It can be concluded, from the given results, that very high compression ratios are achievable by the proposed algorithms. The additional errors introduced by the quantization of wavelet coefficients are very small even when the resolution is reduced to a half (8—bit). In most cases the increase in errors is negligible for the 12-bit resolution, and around 1 to 2 percent for the 8-bit resolution. A significant jump in errors sometimes happens when the compression ratios reach 80 to 100. It should be noted that in almost all the cases — even for the highest compression ratios, the overall relative 61  Chapter 4:  Other Applications of Wavelets in Paper Making Process  error stayed below one, suggesting that even for those cases the achieved accuracy of estimation after compression was better than obtained by the E X P O algorithm without compression. Adaptive methods (MultiMAD) have given much better results (Table 4.20) due to their ability to adjust the thresholds at different levels. The method moved to higher thresholds at some resolution levels so setting more wavelet coefficients to zero (see Table 3.10 for an example).  Case  compress, ratio  8-bit  12-bit  16-bit J compress Jexpo  compress, ratio  Jcompress Jexpo  compress, ratio  Jcompress Jexpo  1.  41.18  0.86  54.91  0.86  82.37  0.87  2.  51.56  1.06  68.74  1.06  103.11  1.07  3.  38.88  .0.91  51.84  0.91  77.77  0.91  4.  51.56  1.06  68.74  1.06  103.11  1.07  5.  36.63  0.89  48.84  0.89  73.26  0.91  6.  53.13  1.07  70.84  1.07  106.26  1.08  7.  13.51  1.11  18.02  1.11  27.03  1.11  Table 4.21 Bump data: summary of compression results  The test data set created by the M V M is used to show the visual effects of the compression. This data set has been chosen because of its richness in M D and C D frequencies yielding relatively smaller compression ratios (Table 4.20). The results described are achieved by the filter #6 (compression ratio of 81.5:1). The energy distribution of the signal is shown in Figure 4.1(a), and is very close to normal distribution. Figures 4.1(b) and 4.1(c) depict, respectively, wavelet coefficients before and k  after the thresholding has been done (Note different scales in (a), (b) and (c)!). It is obvious, from the given diagrams, how wavelet compression operates. The wavelet coefficients are so concentrated in the quantization amplitude levels 50 to 75 that it had to be "zoomed in" to be able to see the low-frequency (high magnitude/energy) wavelet coefficients. The quantization level 62 corresponded to the wavelet coefficient of the zero magnitude. There were 2194 coefficients at this level. Even more impressive are the results after the thresholding was done. A l l the coefficients that used to be grouped around the zero value are now reduced to zero. There is only one peak at the level 62 and 8049 coefficients have that value. This compression resulted in only 143 nonzero coefficients. 62  Chapter 4:  Other Applications of Wavelets in Paper Making Process  Raw data at 8-bit resolution _  80-  i  i  | 60-  *8 20o U 50  100 150 Wavelet coefficients at 8-bit resolution  200  250  100 150 Quantization level Figure 4.1: M V M data set: energy distribution  A n alternative way to show the ability of wavelets to compress a given 2D signal can be seen in Figure 4.2. This figure depicts the wavelet coefficient distribution in a two-dimensional image. Each small black square represents one nonzero wavelet coefficient. The meaning of different regions, marked in Figure 4.2, can be understood by comparison to the diagram given in Figure 4.3 which is related to the calculation of the two-dimensional D W T (see Figure 2.3). The entire rectangle would be occupied by nonzero coefficients prior to compression.  63  Chapter  4:  Other Applications of Wavelets in Paper Making Process  Location of non-zero wavelet coefficients 60 50 40 30 20 10  20  60 Sensors  40  80  100  120  Figure 4.2: MVM.data set: distribution of nonzero wavelet coefficients.  60 50  f  r  1  1  LH  HH  40 30  f  2  20  f LH f 3  10 -  3  LL  f  LH  HH  f  f  3  f  4  HH  3  , HL 20  40  60  i 80  Figure 4.3: 2D DWT: coefficient storage.  64  1  HL i 100  i 120  Chapter  4:  Other Applications of Wavelets in Paper Making Process  The two-dimensional wavelet transform decomposes a signal into a low-frequency (average) signal (fix) and three high-frequency (detail) signals which are directionaly sensitive [18]: fLH — emphasizes horizontal, CD-position dependant, features, fnL — emphasizes vertical, time-dependant, features, fHH — emphasizes diagonal, mutually-dependant, features. This directional sensitivity is an artifact of the frequency ranges they contain. It can be seen (Figure 4.2) that only a few nonzero coefficients lay outside of the f \ L frequency component. It can be confirmed, from Figures 4.4 and 4.5, that the errors are indeed very small and visually they appear as a kind of fuzzyness surrounding the estimated profiles.  65  Chapter  4:  Other Applications of Wavelets in Paper Making Process  Chapter 4:  SUE3S  Other Applications of Wavelets in Paper Making Process  sireos  Chapter  4:  Other Applications of Wavelets in Paper Making Process  Compression has been also tested on the industrial data set given in Section 3.3.2.  The  compression ratios and the errors are given in Table 4.22. The filters are given in the same order as in Table 3.15. In this case the errors (Eq. 3.3) were calculated based on the estimated signal obtained by the same filter since the signal is not known. The errors were very small — even the largest one (case #1) corresponds to an average per sample error of 12.42/(60 * 118) = 0.0018 in comparison to the standard deviation of the estimated signal which was around 0.7. Figures 4.6 to 4.9 are obtained by the filter #3 from Tables 3.15 and 4.22 at the 8—bit resolution. It is impossible to see any difference between compressed and non-compressed images. There is only one wavelet coefficient (Figure 4.7) left outside of the 111 region. Based on the results shown, it is easy to see that the proposed compression method, which is a natural part of the estimation process, yields high compression with almost no additional errors. Different compression strategies may be developed to assure the smallest possible error and the highest compression ratio, based on the use to be made of the stored data. If data is to be used as a simple image-record of the production, good results may be achieved by simply selecting 2—5% of the largest wavelet coefficients before thresholding. This would produce a constant compression ratio. If the record showing the estimated profiles is to be kept, then the data processing described in this section can be used (estimation followed by compression, using the same wavelet coefficients). If there is a need to keep data with high-frequency components still present, thresholding may be done by selecting the threshold conservatively (i.e. one half of the one automatically selected by the estimator) and by using non-adaptive methods. This filtering would remove a portion of the noise but would still achieve a useful compression ratio.  68  Chapter  4:  Other Applications of Wavelets in Paper Making Process  (a) Raw data at 8-bit resolution  300  (c) Thresholded wavelet coefficients at 8-bit resolution  Quantization level Figure 4.6: Pulp machine moisture data: energy distribution  Case  compress,  compress,  compress, ratio  8-bit  12-bit  16-bit Jcompress  ratio  Jcompress  ratio  Jcompress  1.  11.13  0.06  17.77  0.19  26.66  12.42  2.  13.62  0.01  18.16  0.42  27.24  2.37  3.  13.72  0.00  18.24  0.06  27.44  0.89  Table 4.22 Pulp machine moisture data: summary of compression results 69  Chapter 4:  Other Applications of Wavelets in Paper Making Process  Location of non-zero wavelet coefficients  20  40  60 Sensors  80  100  120  Figure 4.7: Pulp machine moisture data: distribution of nonzero wavelet coefficients.  70  Chapter 4:  (%) ajnjsic-w  Other Applications of Wavelets in Paper Making Process  (%) airnsiow  71  Chapter  4:  72  Other Applications of Wavelets in Paper Making Process  Chapter  4:  Other Applications of Wavelets in Paper Making Process  4.2 Wavelet Estimation for Control Purposes The wavelet estimator presented in Chapter 3 produced very good results with low estimation errors and at high execution speeds. Those estimates, in the form of two- and three-dimensional images can provide a paper machine operator with a much better view of the production quality. The operator acts on the basis of this improved information leading ultimately to better control and paper quality. The estimator, as given in Chapter 3, worked in the batch mode on a relatively large number of data samples (more then 50 scans and more then 110 sensor measurements per scan) and was quite suitable for the above purposes. Some modifications to the proposed algorithms are needed for the estimation to become more appropriate for use in real-time control. Those modifications concentrate on more practical issues such as the speed of M D profile updates and the handling of edge effects. This is mainly a topic for future research. One possible approach to this problem is given as follows. The C D profiles should be estimated as shown in Chapter 3. The number of scans that need to be used for each C D profile calculation depends on the filter length and it should be at least in the range of 2 to 3 filter lengths. The most recent estimates of the signal are the most important because these are the samples needed for the controller to respond at the highest rate. Hence, better handling of those data samples may be needed (i.e. using Meyer's boundary wavelets or Dyadic boundary wavelets). After the C D profile is calculated (once per scan, batch mode) the M D profiles may be estimated by subtracting the C D profile from the raw data and using the result as an input to a one-dimensional wavelet de-noising procedure. This approach may lead to M D estimation at each sample point instead of once per scan. A n example of the M D profile estimation at each sample and its comparison to the scan average method is given in Figure 4.10. The data set used in this example was given in Section 3.2.2 (Smooth Surfaces). The true M D profile, its scan average estimate, as well as the high resolution estimates achieved by the proposed method are shown in Figure 4.10(b). The raw M D profile is shown in Figure 4.10(a) for comparison. The same method was applied to the pulp machine data (Section 3.3.2). The results as shown in Figure 4.11(a)-(b). The wavelet estimator will thus enable more frequent M D estimates. The M D controller will then be able to operate with reduced time delay and possibly with increased bandwidth. A second possible benefit is in C D control, where the profile estimate, after decomposition into wavelet components, can be matched to C D actuator spacing and response characteristic. 73  Chapter 4:  Other Applications of Wavelets in Paper Making Process  (a) Raw MD profile  (b) MD profile comparason  True MD ScanAvg.  -lh -  ID WAV 10  40  30 Scans  20  50  60  Figure 4.10: Smooth surfaces: high resolution MD estimation (a) Raw MD profile  (b) MD profile comparason 10  • ScanAvg. 1D_WAV 10  20  30 Scans  40  50  60  Figure 4.11: Pulp machine moisture data: high resolution MD estimation  74  Chapter  5:  Conclusion  Chapter 5 Conclusion The thesis describes a new approach to paper machine process data analysis using onedimensional and two-dimensional discrete wavelet transforms. These techniques have been adapted from a general theory that has been developed in recent years on the application of wavelets to signal analysis. Application areas in which the theory was first applied have included image processing and bandwidth compression for communications. Two main applications of the discrete wavelet transform have been analyzed in this thesis. First, an analysis of the use of wavelets for processing scanned data representing basis weight and moisture variations on a paper machine has been carried out. It has been shown that wavelets are effective for the detection of process signals in noisy data, so leading to better estimation and visualization of the machine direction and cross machine variations in process data. The second main application of the method has been to allow significant compression of the process data without diminishing the ability to reconstruct accurate profiles. It has been shown that the compression method can be embedded into the estimation algorithm, producing excellent results without a major expense in the computation time. It has been shown that, in both applications, the new methods produce results superior to the industrially accepted procedures.  For appropriate choice of wavelets, profile estimates are  improved over those obtained using exponential filtering or other standard analysis methods. The data compression technique presents a new concept in paper machine data analysis and the author is not aware of any previous references to this subject. The ability to reduce data storage requirements is of importance in mill-wide process monitoring systems. A comprehensive analysis of the proposed algorithms has been carried out on a variety of simulated data sets for which the true process variations are known. Industrial data has also been analyzed and it is apparent that the method had many desirable characteristics. The following properties have been shown for the wavelet estimator: •  The estimator reduces the error of estimation compared to the existing algorithms, when applied to test data. 75  Chapter  •  5:  Conclusion  When applied to industrial data, the method presents an excellent visual interface to the paper making process. Operator data presentations are an important part of process control systems, and it is believed that wavelet analysis gives excellent visual insight into process variations, insight significantly better than is available in current systems.  •  The wavelet decomposition can be closely related to cross machine actuator responses. It is expected that a cross machine control method based on an appropriate wavelet analysis will allow profile upsets to be matched to actuator settings.  •  The algorithm is computationally fast, which allows real time analysis. The need for a recursive form may be avoided. The method has been shown to be applicable to basis weight and moisture estimation using data from a scanning sensor.  It is to be expected that it could be used for the  analysis of other process variables. •  Unlike many estimation methods based on least squares techniques, wavelet analysis is robust in its performance and does not require additional tuning once it has been properly installed.  The wavelet compression algorithm has been shown to have the following advantages: •  The method is very effective with compression ratios greater than 25:1 (up to 100:1) and with the error of estimation, after decompression, remaining lower than the corresponding error for an exponential filter based algorithm without compression.  •  Data compression is an integral part of the estimator procedure, and therefore requires very little additional computation time.  •  The algorithm provides an effective tool for the visualization of production variations over large time periods. It follows that it can be used to search for the possible causes of disturbances and faults.  •  The algorithm can be adjusted to provide a minimum error or a constant compression ratio, allowing it to be tailored to particular applications.  •  Data compression will become increasingly important with the introduction of new, high speed sensors that generate large volumes of data. A related advantage is the ability to reduce time for data transfers over the local area networks used in modern machines.  76  Chapter  5:  Conclusion  There are many possible extensions to the research presented here. Steps that might be taken to improve the effectiveness of the algorithms are described in the following: •  There is a need for better handling of the edge effects for the given D W T algorithms.  •  For real time data analysis it is desirable that recursive versions of the algorithms replace the batch procedures.  •  A n organized procedure for incorporating knowledge of the process and noise models into the estimation procedures is needed. As the theory associated with wavelet decomposition advances, there is a need to evaluate new approaches. For example, there is a case for exploring new wavelet techniques such as wavelet packet analysis [38] and second generation wavelets [35].  •  It is likely the use of wavelet estimation algorithms in automatic control applications will allow sensor and actuator spatial and temporal bandwidths to be matched.  •  The existing threshold estimators that were used and play a critical role in data compression and noise elimination can still be improved. In some cases (for example the SURE technique), the results were not as expected and there is a need for further development of this part of the method.  77  Chapter  5:  Conclusion  References  [1]  J. Buckheit, S. Chen, D . Donoho, I. Johnstone, and J. Scargle. Wavelab 0.655, a library of matlab routines for wavelet analysis, ftp://playfair.stanford.edu/pub/wavelab, 1995.  [2]  S.-C. Chen. Kalman filtering applied to sheet measurement. In 7th American Control Conference, pages 643-647, Atlanta, Georgia, 1988.  [3]  A . Cohen, I. Daubechies, and J. C. Feauveau. Biorthogonal bases of compactly supported wavelets. Comm. Pure Applied Math., 1992.  [4]  E. Dahlin. Computational methods of a dedicated computer system for measurement and control on paper machines. In 24th Engineering Conference, TAPPI, pages 62.1—62.42, San Francisco, California, Sep 1969.  [5]  I. Daubechies. Orthonormal bases of compactly supported wavelets. Comm. Pure Applied Math., XLI(41):909-996, Nov. 1988.  [6]  I. Daubechies. Orthonormal bases of wavelets with finite support - connection with discrete filters. In J. M . Combes, A . Grossman, and P. Tchamitchian, editors, Wavelets, TimeFrequency Methods and Phase Space, pages 38-66. Springer-Verlag, Berlin, 1989. Proceedings of International Colloquium on Wavelets and Applications, Marseille, France, December 1987.  [7]  I. Daubechies. The wavelet transform, time-frequency localization and signal analysis, Sept. 1990.  [8]  I. Daubechies. Ten Lectures on Wavelets. SIAM, Philadelphia, PA, 1992. Notes from the 1990 C B M S - N S F Conference on Wavelets and Applications at Lowell, M A .  [9]  I. Daubechies. Orthonormal bases of compactly supported wavelets II. Variation on a theme. SIAM Journal of Mathematical Analysis, 24(2):499-519, March 1993.  [10]  D . L . Donoho. Nonlinear wavelet methods for recovery of signals, densities, and spectra from indirect and noisy data. In Proceedings of Symposia in Applied Mathematics, volume 00, pages 173—205. American Mathematical Society, 1993. 78  Chapter  [11]  5:  Conclusion  D. L . Donoho. De-noising by soft-thresholding. IEEE Trans. Inform. Theory, 41(3):613—627, May 1995.  [12]  D . L . Donoho and I. M . Johnstone. Minimax estimation via wavelet shrinkage. Technical Report 402, Stanford University, Department of Statistics, July 1992.  [13]  D . L . Donoho and I. M . Johnstone. Adaptating to unknown smoothness via wavelet shrinkage. J. Am. Stat. Ass., 1993. To appear, Also Tech. Report 425, Department of Statistics, Stanford University, June, 1993.  [14]  D . L . Donoho and I. M . Johnstone. Ideal spatial adaptation via wavelet shrinkage. Biometrika, 81:425-455, 1994. Also Tech. Report 400, Department of Statistics, Stanford University, July, 1992.  [15]  G . A . Dumont, M . S. Davies, K . Natarajan, C. Lindeborg, F. Ordubadi, Y . Fu, K . Kristinsson, and A . Jonsson. A n improved algorithm for estimating paper machine moisture profiles using scanned data. In 30th IEEE Conference on Decision and Control, Brighton, England, Dec 1991.  [16]  R. A . Gopinath, H . Guo, R. Hindman, M . Lang, J. E . Odegard, and D . Wei. Rice university wavelet toolbox, version 2.01. http://www-dsp.rice.edu/software/RWT, Apr 1994.  [17]  A . Haar. Zur theorie der orthogonalen funktionensysteme. Mathematische Annalen, 69:331— 371, 1910. Also in PhD thesis.  [18]  M . Hilton, B . Jawerth, and A . Sengupta. Compressing still and moving images with wavelets. Multimedia Systems, 3(2), 1994. To appear.  [19]  A . K . Jain. Fundamentals of digital image processing. Prentice Hall, Englewood Cliffs, N J , 1989.  [20]  B . Jawerth and W. Sweldens. A n overview of wavelet based multiresolution analyses. Technical Report 1993:1, Industrial Mathematics Initiative, Department of Mathematics, University of South Carolina, 1993. Submitted to SLAM Review.  [21]  I. M . Johnstone and B . W. Silverman. Wavelet threshold estimators for data with correlated noise. Technical report, University of Bristol, U K , Statistics Department, Sept. 1994.  [22]  I. M . Jonsson. Estimation and identification of moisture content in paper. Master's thesis, The University of British Columbia, Aug. 1991.  Chapter  [23]  5:  Conclusion  M . Lang, H . Guo, J. E . Odegard, C. S; Burrus, and R. O. Wells, Jr. Noise reduction using an undecimated discrete wavelet transform. In To appear in IEEE Signal Processing Letters, 1995.  [24]  P.-G. Lemarie and Y . Meyer. Ondelettes et bases hilbertiennes. Rev. Mat.Tberoamericana, 2:1-18, 1986.  [25]  C . Lindeborg. A nonlinear algorithm for the estimation of moisture variation characteristics in the paper process. In Proc. of 17th International Conference BIAS-81, volume 3, pages 139-158, Milan, Italy, Oct 1981.  [26]  S. Mallat. A theory for multiresolution signal decomposition: The wavelet representation. IEEE Trans, on Pattern and Machine Intelligence, 11(7):674-693, 1989.  [27]  D . Marr. Vision: A Computational Investigation into the Human Representation and Processing of Visual Information. W. H . Freeman, 1982.  [28]  Y . Meyer. Ondelettes, functions splines et analyses graduees. Lectures at Univ. of Torino, Italy, 1986.  [29]  Y . Meyer. Ondelettes. filtres miroirs en quadrature et traitement numerique de l'image. Technical report, Cahiers de Mathematiques de la Decision, Paris, France, 1989. Text from a conference on 22 February 1989 at the University of Autonooma of Madrid.  [30]  S. T. Morgan. Estimation and identification for machine direction control of basis weight and moisture. Master's thesis, The University of British Columbia, June 1994.  [31]  G. P. Nason. Wavelet regression by cross-validation. Technical Report 447, Stanford University, Statistics Department, Mar. 1994.  [32]  K . Natarajan, G . A . Dumont, and M . S. Davies. A n algorithm for estimating cross and machine direction moisture profiles for paper machines. In IFFAC/FORS Symposium, pages 27—31, Beijing, P R C , Aug 1988.  [33]  W. K . Pratt. Digital Image processing. Wiley, New York, 1978.  [34]  O. Rioul and M . Vetterli. Wavelet and signal processing. IEEE Signal Processing Magazine, 8(4):14-38, Oct. 1991. 80  [35]  W. Sweldens. The lifting scheme: a construction of second generation wavelets. Technical Report 1995:6, Industrial Mathematics Initiative, Department of Mathematics, University of South Carolina, 1995.  [36]  B . F. Taylor. Optimum separation of machine-direction and cross-direction product variations. Tappi Journal, pages 87-92, Feb 1994.  [37]  X . G . Wang, G. A . Dumont, and M . S. Davies. Modelling and identification of basis weight variations in paper machines. IEEE Transactions on Control System Technologies, June 1993.  [38]  M . V . Wickerhauser. Smooth localized orthonormal bases. CRAS, 315, 1993.  81  Appendix A Exponential Multiple Scan Trending (EXPO) Algorithm The equations describing the algorithm (Dahlin [1]) are:  (A.1) s (m) = p{D (m) - Y(m)) + (1 - p)s (m - 1) n  n  n  where: D (m) n  is the measured value (moisture content or basis weight) at scan m and C D position n;  Y(m) is the mean value o f scan m; N is the total number of samples (data boxes) in a C D ; s (m) n  is the filtered C D profile of scan m at data box n; and  p is the exponential filter weighting factor, 0 < p < 1; The weighting factor is usually chosen from the region p < 0.3. This eliminates the M D variations from the C D profiles but also slows down the C D profile estimation.  82  Appendix B Estimation and Identification of Basis Weight and Moisture Content (EIBMC) Algorithm B.l Process Models Process model used to describe the moisture variation in the ELBMC algorithm is based on Lindeborg's [3] moisture variation model (MVM). It has been shown [6] that the same model, with some modifications, can be applied to basis weight. In this thesis only the necessary equations, describing the process models and the ELBMC algorithm are given, For the detailed discussion one can refer to the original papers [3, 6, 2, 4].  Machine direction wet-end "  reel  Cross direction  k+1 k-1  Figure B . l : Data indexing  B.l.l Moisture The M V M for dry-end variations is Vk=Pk + £ + BPk>k + v  k  (B.l)  where n is data box number (CD-position); k is sample number or time instant (kT); yl is the difference between measured value at C D position n and time k, and the variable set point; B is a measure of the non-linearity between the C D and M D variations; 83  pi is pseudo-profile at position n and time k; ut is M D value at time k; vk ~ N(0, R) is the measurement noise and model mismatch; It is assumed that P  k  =  {p\,v\, • • • ,Pk]  a  n  d  change slowly [3]. The first order model of  B  the M D variations is given [5] as  = u +&  Uk  (B.2)  where u is the mean moisture content in the M D ; £jk is a zero mean stochastic process given by:  & = a£k-i +Wk,  (B.3)  Wk ~ N(0, q)  The parameter a and the variance q depend on the paper machine. The model in state space form is given by  x i = Ax + W k+  k  k  (B.4) Vk=Pk+ C%x + v k  where  "1 0"  11 '  A =  Xk  h. B.1.2  k  0 (B.5)  0 a.  Cl = [l  \\{l + B l) V  Basis Weight  It has been shown [6] that the previously described moisture model ( B . l . l ) can be modified to accommodate basis weight estimation by increasing the order of the disturbance (second order A R M A ) , setting B=0 and absorbing the state u into the A R M A model. The A R M A model is given as  ( l - aiq'  1  - a q~ )u 2  2  = ( l + biq' 4- hq~ )w 1  k  2  k  (B.6)  Thus, the state space model for the basis weight becomes:  x +i = Ax + W k  k  Uk  = Cxk +  Vk=Pk+Uk  84  k  (B.7)  Wk  + v  k  where  A  a\  1  a  0  Wk 2  (B.8)  C = [l 0] B.2 The E I B M C Algorithm  B.2.1 Moisture The moisture version of the ELBMC algorithm is  k  Uk\k-i = I  k  k  = 1+  n  (B.ll)  Bu \k-l k  Vk\k-i = Pk-i + V " = -V  (B.10)  l]^fc|fc-i  1  ip  (B.9)  Ax _i\ _i  Xk\k-i =  + p-s(v y  -  n  av r (y -y ^) n  n  k  (B.13)  n  k  klk  Vl = p g - i +  V»  i + = ( l + B p j ) [1  1]  Sfc-i(C^) "  (B..12)  Cfc-i*fc|fc-i  T  •q?E _ (q?) -rJJ  (B.14) (B.15) (B.16)  T  f c  1  (B.17) (B.18)  Wfc|fc = [1  (B.19)  1 ]x \ k  k  where Q = E{W W ] r  k  k  =  [91 0 .0  <?J  The equations are given in the order of their execution. 85  (B.20)  B.2.2 Basis Weight The basis weight version of the E I B M C algorithm is  Xk\k-i =  (B.21)  Ax -\\k-i k  (B.22)  u \ -\=Cx \ _ k  k  k  k  1  03.23)  Vk\k-1 = Vk-\ + Uk\k-\ V  n  =  \yn  _  Q  ( ") y  0 _  2 +  fitynf  1+V  A  .  n  •'  (B.25)  v  n  1+ Vk — Pk-\ + Ki CT _iC  (B.24)  _  (B.26) T  f c  + .R (B.27)  Xk\k = Xk\k-i + K*k {Vk ~ {Pk + Uk\k  T  k  k  k  =  b\  hb  hh  b\  2  The equations are given in the order of their execution.  86  k  (B.28) (B.29)  = Cxk tk\k  where Q = E{W W )  Cx \ _i))  (B.30)  References [1]  E . Dahlin. Computational methods of a dedicated computer system for measurement and control on paper machines. In 24th Engineering Conference, TAPPI, pages 62.1—62.42, San Francisco, California, Sep 1969.  [2]  G . A . Dumont, M . S. Davies, K . Natarajan, C. Lindeborg, F. Ordubadi, Y . Fu, K . Kristinsson, and A . Jonsson. A n improved algorithm for estimating paper machine moisture profiles using scanned data. In 30th IEEE Conference on Decision and Control, Brighton, England, Dec 1991.  [3]  C. Lindeborg. A process model of moisture variations. Pulp and Paper Canada, 4:T142—147, 1993.  [4]  S. T. Morgan. Estimation and identification for machine direction control of basis weight and moisture. Master's thesis, The University of British Columbia, June 1994.  [5]  K . Natarajan, G . A . Dumont, and M . S. Davies. A n algorithm for estimating cross and machine direction moisture profiles for paper machines. In IFFAC/FORS Symposium, pages 27—31, Beijing, P R C , Aug 1988.  [6]  X . G . Wang, G . A . Dumont, and M . S. Davies. Modelling and identification of basis weight variations in paper machines. IEEE Transactions on Control System Technologies, June 1993.  87  Appendix C  Wavelet Filter Coefficients  Haar 0.70710678118655 0.70710678118655 Daubechies(4) 0.48296291314453 0.83651630373781 0.22414386804201 -0.12940952255126 Daubechies(6) 0.33267055295008 0.80689150931109 0.45987750211849 -0.13501102001025 -0.08544127388203 0.03522629188571 Symmlet(6) 0 . 01540410932734 0.00349071208433 -0.11799011114841 -0.04831174258600 0.49105594192764 0 .78764114102879 0 . 33792942172824 -0.07263752278660 -0.02106029251270 0.04472490177075 0. 00176771186440 -0.00780070832477 Symmlet(8) 0.00188995033290 -0.000302.92051455 -0.01495225833679 0.00380875201406 0.04913717967348 -0.02721902991681 -0.05194583810788 0 .36444189483596 0 .77718575169975 0.48135965125920 -0.06127335906791 -0.14329423835105 0.00760748732528 0.03169508781035 -0.00054213233164 -0.00338241595136  Coiflet(2) 0.01638733646360 -0.04146493678196 -0.06737255472228 0.38611006682299 0.81272363544940 0.41700518442367 -0.07648859907867 -0.05943441864674 0.02368017194645 0. 00561143481942 -0.00182320887071 -0.00072054944537 Coiflet(3) -0.00379351286449 0.00778259642733 0.02345269614184 -0.06577191128186 -0.06112339000267 0.40517690240961 0.79377722262562 0.42848347637762 -0.07179982161931 -0.08230192710689 0.03455502757306 0.01588054486362 -0.00900797613666 -0.00257451768875 0.00111751877089 0.00046621696011 -0.00007098330314 -0.00003459977284  

Cite

Citation Scheme:

    

Usage Statistics

Country Views Downloads
United States 15 4
China 12 34
India 12 0
Russia 6 0
Canada 5 0
Japan 3 0
Libya 2 0
Malaysia 2 0
France 2 0
Sweden 1 0
Egypt 1 0
Germany 1 0
City Views Downloads
Unknown 18 0
Beijing 8 0
Patel 5 0
Shenzhen 4 34
Ashburn 4 0
Tokyo 3 0
Ridgeland 3 0
Alangulam 2 0
Potsdam 2 0
Ottawa 2 0
Dundalk 2 0
Tripoli 2 0
Stockholm 1 0

{[{ mDataHeader[type] }]} {[{ month[type] }]} {[{ tData[type] }]}
Download Stats

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0065195/manifest

Comment

Related Items