UFDC Home  myUFDC Home  Help 



Full Text  
PAGE 1 1 AUTOMATED STRAIN ANALYSIS SYSTEM: DEVELOPMENT AND APPLICATIONS By WEIQI YIN A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 2009 PAGE 2 2 2009 Weiqi Yin PAGE 3 3 To my pare nts Guoyin and Jine, and my wife, Min PAGE 4 4 ACKNOWLEDGMENTS I sincerely thank my advisor, Dr. Peter Ifju, for his support, guidance, knowledge and friendship. It is very fortunate for me to work with him. A lso, I sincerely thank my committee members, Dr. Bhavani Sankar, Dr. Raphael Haftka, Dr. Youping Chen and Dr. Mang Tia for participating and evaluating my research work. I would like to e specially thank Tzuchau Chen, my colleague in the Experimental Stress Anal ysis lab for his continuous help o n experi ments and shared knowledge. I would also like to thank my other lab colleagues, Nancy Strickland, Vijay Jagdale and Mulugeta Haile, for their shared knowledge and help. I am very grateful to Christian Gogu for his shared knowledge on material property identification using full field measurement. PAGE 5 5 TABLE OF CONTENTS ACKNOWLEDGMENTS .................................................................................................. 4 page LIST OF TABLES ............................................................................................................ 8 LIST OF FIGURES .......................................................................................................... 9 ABSTRACT ................................................................................................................... 14 CHAPTER 1 INTRODUCTION .................................................................................................... 16 Background ............................................................................................................. 16 Literature Review .................................................................................................... 20 Intensity Based Analysis .................................................................................. 21 Phase Extraction Methods ................................................................................ 22 Phase Unwrapping ........................................................................................... 23 Strain Calculation ............................................................................................. 25 Objectives ............................................................................................................... 26 2 APPROACH ............................................................................................................ 28 Pre Processing of Fringe Pattern ............................................................................ 28 Image Processing in the Spatial Domain .......................................................... 29 Smoothing linear filters .............................................................................. 31 Order statistics filters ................................................................................. 31 Image Processing in the Frequency Domain .................................................... 32 Ideal low pass filters ................................................................................... 32 Gaussian low pas s filters (GLPF) ............................................................... 33 Summary .......................................................................................................... 33 Phase Shifting ......................................................................................................... 34 Theory of Phase Shifting .................................................................................. 35 Experimental Setup .......................................................................................... 37 Stage for Phase Shifting ................................................................................... 38 E limination of Miscalibration (N+1 Phase Shifting) ........................................... 40 Phase Unwrapping ........................................................................................... 40 Itohs method ............................................................................................. 41 Quality guided phase unwrapping .............................................................. 44 Improved quality guided phase unwrapping ............................................... 49 Repair of inconsistent areas ....................................................................... 52 Detection of inconsistent areas .................................................................. 54 Summary .......................................................................................................... 57 Displacem ent and Strain Calculation ...................................................................... 57 Traditional Method ............................................................................................ 57 PAGE 6 6 Global Surface Fit Based Strain Calculation ..................................................... 61 Local Surface Fit Based Strain Calculation ...................................................... 65 Summary .......................................................................................................... 68 Optical/Digital Fringe Multiplicatio n (O/DFM) Method ............................................. 69 Theory of O/DFM Method ................................................................................. 69 Fringe Skeleton Detection ................................................................................ 71 Fringe Thinning ................................................................................................ 74 Fringe Order Assignment ................................................................................. 75 Fringe Order Interpolation ................................................................................ 76 Strain Calculation ............................................................................................. 77 Summary .......................................................................................................... 78 Summary ................................................................................................................ 78 3 AUTOMATED STRAIN ANALYSIS SYSTEM ......................................................... 80 Platform .................................................................................................................. 80 Functions of the System ......................................................................................... 80 Demonstration of the System .................................................................................. 81 Validation ................................................................................................................ 88 4 point Bending Test ........................................................................................ 88 Openhole Tension Test ................................................................................... 91 Summary ................................................................................................................ 94 4 APPLICATIONS OF AUTOMATED STRAIN ANALYSIS SYSTEM ....................... 96 Residual Stress Characterization of Plain Woven Composites ............................... 96 Introduction ....................................................................................................... 96 C ure Reference Method ................................................................................... 99 Result of Residual Strains from Experiment ................................................... 101 Finite Element Model Description ................................................................... 103 Result of Residual Strain and Stress from FE Model ..................................... 107 Summary ........................................................................................................ 110 Shrinkage Measurement of Concrete Material ...................................................... 111 Introduction ..................................................................................................... 111 Specimen Preparation .................................................................................... 112 Result for Cement without Aggregate ............................................................. 113 Result for Cement with Coarse Aggregate ..................................................... 119 Summary ........................................................................................................ 121 Material Property Identification of Laminate Using Full Field Measurements ........ 121 Introduction ..................................................................................................... 121 Experiment ..................................................................................................... 123 Fringe Patterns and Result ............................................................................. 125 Summary ........................................................................................................ 127 Summary .............................................................................................................. 127 5 CONCLUSIONS AND SUGGESTED FUTURE WORK ........................................ 128 PAGE 7 7 Conclusions .......................................................................................................... 128 Suggested Future Work ........................................................................................ 129 APPENDIX A BILINEAR FUNCTION FOR LOCAL SURFACE FIT BASED STRAIN CALCULATION ..................................................................................................... 130 B CUBIC FUNCTION FOR LOCAL SURFACE FIT BASED STRAIN CALCULATION ..................................................................................................... 132 LIST OF REFERENCES ............................................................................................. 134 BIOGRAPHICAL SKETCH .......................................................................................... 144 PAGE 8 8 LIST OF TABLES Table page 2 1 Procedure of Itohs method of phase unwrapping .............................................. 42 2 2 Procedure of quality guided phase unwrapping .................................................. 48 2 3 CPU time vs. image size .................................................................................... 50 2 4 Improved procedure of quality guided phase unwrapping .................................. 51 2 5 CPU time vs. division size .................................................................................. 52 2 6 Procedure of detection and repair of inconsistent areas ..................................... 56 3 1 List of functions in the software .......................................................................... 80 3 2 Concentration for the openhole aluminum test .................................................. 94 4 1 Material properties used in FE model for plain weave composite RVE ............ 106 4 2 Manufacturers specifications and properties found by Noh (2004) based on a four points bending test .................................................................................... 123 4 3 Mean values and coefficient of variation of the identified posterior distribution from the Moir interferometry experiment ......................................................... 126 PAGE 9 9 LIST OF FIGURES Figure page 1 1 Four beam Moir Interferometer schematic ....................................................... 17 1 2 Four beam Moir Interferometer ......................................................................... 17 1 3 Typical fringe patterns (scribed circle is 1inch diameter) ................................... 18 2 1 Illustration of large random noise existing in fringe pattern. A) Original fringe pattern. B) Intensity dis tribution along the lines ................................................. 28 2 2 Illustration of large random noise existing in wrapped phase. A) Wrapped phase. B) Wrapped phase distribution along the lines ....................................... 29 2 3 Spatial filtering. The magnified drawing shows a 33 mask and the image section directly under it ....................................................................................... 30 2 4 Wrapped phase map from 4 fringe patterns filter ed via 55 averaging filter ....... 34 2 5 Wrapped phase map from 4 fringe patterns filtered via GLPF filter .................... 34 2 6 Generalized twobeam interferometer ................................................................ 35 2 7 Phasor diagram .................................................................................................. 35 2 8 Setup for experiment .......................................................................................... 38 2 9 Stage designed for phase shifting ...................................................................... 39 2 10 Calibration of the stage ....................................................................................... 39 2 11 Wrapped phase and unwrapped phase .............................................................. 4 0 2 12 Wrapped phase and unwrapped phase .............................................................. 43 2 13 Wrapped phase and unwrapped phase A) Wrapped phase B) Unwrapped phase obtained using Itohs method ................................................................... 44 2 14 Integration direction of residue calculation ......................................................... 45 2 15 Residues detected .............................................................................................. 46 2 16 Wrapped phase map and its quality map. A) Wrapped phase map. B) PDV quality map ......................................................................................................... 47 PAGE 10 10 2 17 Process of Quality guided phase unwrapping A) 100 pixels unwrapped. B) 1000 pixels unwrapped. C ) 10000 pixels unwrapped D ) 50000 pixels unwrapped. E ) 80000 pixels unwrapped. F ) Final unwrapped phase ............... 48 2 18 Unwrapped phase map. A) by Itohs method B) by quality guided phase unwrapping ......................................................................................................... 49 2 19 Unwrapped phase map. A) by Itohs method B) by quality guided phase unwrapping ......................................................................................................... 51 2 20 Removal of false fringe manually. A) Wrapped phase with false fringe. B) Zoomed in the area with false fringe. C) Unwrapped phase with big error around the false fringe. D) Wrapped phase after removal of false fringe. E) Unwrapped phase with false fringe removed. ..................................................... 52 2 21 Detection and repair of inconsistent areas. A) Original unwrapped phase. B) Quality map. C) Inconsistent areas detected. D) Extended inconsistent areas. E) Repaired unwrapped phase map. F) Difference of unwrapped phase before and after repair ............................................................................. 54 2 22 Detection of inconsistent areas ........................................................................... 55 2 23 4 point bending. A) Test s cheme B) 4 consecutive phase shifted fringe patterns .............................................................................................................. 59 2 24 Calculation of normal strain xx A) Gage length = 10 pixels B) Gag e length = 30 pixels C) Strain along y=10 pixels ............................................................ 60 2 25 Normal strain xx obtained via global surface fit ................................................ 64 2 26 Normal strain xx obtained via local surface fit A) Quadric B) Bilinear ........... 67 2 27 Plot of strain along y=10 pixels ........................................................................... 69 2 28 Explanation of O/DFM A) 2 complementary fringes B) Subtraction of each other C ) Absolute value of B) D ) After binarizing ............................................ 70 2 29 Fringe skeleton det ected using binarizing method for 02II ......................... 72 2 30 Fringe skeleton detected local minimum ............................................................ 72 2 31 Fringe skel eton detected using local minimum detection method for 02II .. 73 2 32 Fringe skeleton after fringe thinning for 02II ............................................... 74 2 33 Fringe order assignment ..................................................................................... 75 PAGE 11 11 2 34 Fringe order interpolation ................................................................................... 76 2 35 Strain calculated. A) By phase shifting. B) By O/DFM ....................................... 77 3 1 Openhole test specimen and fringe patterns ..................................................... 82 3 2 Open image files ................................................................................................. 82 3 3 Coordinate system translating parameters ......................................................... 83 3 4 Noise removal ..................................................................................................... 83 3 5 Define the boundary ........................................................................................... 84 3 6 Quality guided phase unwrapping. A) Unwrapping progress B) Final unwrapped phase ............................................................................................... 85 3 7 Phase smoothing ................................................................................................ 86 3 8 Data extraction. A) Define the data path. B) Data along the data path .............. 87 3 9 Fringe pattern. A) U field. B) V field ................................................................... 88 3 10 Strain maps from automated strain analysis system. A) Strain xx B) Strain yy C) Strain xy ........................................................................................ 89 3 11 Strain maps from FE. A) Strain xx B) Strain yy C) Strain xy ........................... 90 3 12 S train result from manual calculation and analysis system. A) Strain xx from manual calculation. B) Strain xx from the analysis system ................................. 91 3 13 S train map ( yy ) from automated strain analysis system and FE simulation. A) From automated strain analysis system. B) From FE simulation ................... 92 3 14 Strain result from man ual calculation and analysis system. A) Strain xx from manual calculation. B) Strain xx from the analysis system ................................. 93 3 15 Strains along the narr owest section .................................................................... 94 4 1 Representative volume element of the plain woven textile ................................. 98 4 2 Cure reference method ....................................................................................... 99 4 3 4 CRM gratings ................................................................................................. 100 4 4 Fringe patterns. A) U field B) V field ................................................................ 101 PAGE 12 12 4 5 Displac ement field. A) U field B) V field ........................................................... 101 4 6 Strain map. A) Strainxx B) Strain yy ............................................................... 102 4 7 Strain map within one RVE. A) Strain xx B) Strain yy ..................................... 103 4 8 Geometry shape ............................................................................................... 104 4 9 Geometry model of fiber and matrix. A) Fiber B) Matrix .................................. 105 4 10 Strain map from FE modeling. A) Strain xx B) Strain yy ................................. 108 4 11 Residual stress maps from FE modeling. A) Stressxx B) Stressyy .............. 109 4 12 U f ield under 50% RH and 23oC ....................................................................... 114 4 13 V field under 50% RH and 23oC ....................................................................... 114 4 14 Displacement field in the Rectangular Coordinate syst em for day 2 ................. 115 4 15 Displacement field in the Polar Coordinate system for day 2 ........................... 115 4 16 Strain field in the Rectangular Coordinate system for day 2 ............................. 116 4 17 Strain field in the Polar Coordinate system for day 2 ........................................ 117 4 18 Shrinkage from center to expose s urface ......................................................... 118 4 19 Temperature effect ........................................................................................... 118 4 20 Relative humidity effect .................................................................................... 119 4 21 The side and top views of gravel test ............................................................... 119 4 22 Day 5 moir fringe patterns for the gravel test .................................................. 120 4 23 The normal strain analysis for gravel test ......................................................... 120 4 24 Specimen geometry. The specimen material is graphite/epoxy and the stacking sequence 45,45,0s ......................................................................... 123 4 25 Test specimen with grating (1200 lines/mm) .................................................... 124 4 26 Experimental setup ........................................................................................... 124 4 27 Fringe patterns. A) U field B) V field ................................................................ 125 4 28 Displacement maps. A) U field B) V field ......................................................... 126 PAGE 13 13 4 29 Strain maps. A) Strain xx B) Strain yy ............................................................. 127 PAGE 14 14 Abstract of Dissertation Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy AUTOMATED STRAIN ANALYSIS SYSTEM: DEVELOPMENT AND APPLICATIONS By Weiqi Yin December 2009 Chair: Peter G. Ifju Major: Mechanical Engineering I nterferometry has been widely applied for industrial and scientific studies because of its ca pability of fullfield displacement measurement, high sensitivity and high spatial resolution. T he output of m oir is represented as interference fringe patterns which need to be analyzed to obtain the desired physical parameters such as displacement s strain s, etc Numerous algorithms have been developed for fringe pattern analysis however, most of them are focused on specific aspect and hard to be automated, and the strain result is usually a qualitative analysis. I n t his dissertation, the existing techniques were investigated, the appropriate algorithms were adopted, s ome of the selected algorithms were improved, and new algorithms were devel oped. A technique based on phase extraction is investigated, improved and applied for fringe pattern analysis P hase s hifting method uses a series of phase shifted fringe patterns and involves procedures of noise filtering, wrapped phase calculation, phase unwrapping and calculation of displacement and strain. Appropriate digital image processing techniques are applied for noise filtering. P hase shifting algorithm is used to calculate wrapped phase. A quality guided phase unwrapping method is adopted and PAGE 15 15 improved. Global and local surface fit based smoothing algorithms are developed to calculate displacement s and strains. Additionally, a hybrid method based on both intensity and phase information, optical/digital fringe multiplication (O/DFM), is also investigated and developed. O/DFM also uses a series of phase shifted fringe patterns and involves the steps of noise filtering, fringe multiplication, fringe centerline detection, fringe thinning, fringe order assignment, fringe order interpolation, and calculation of displacement and strain. An automated strain analysis system based on the above techniques is developed as a Windows GUI based software in th e Experimental Stress Analysis (ESA) Lab. B y the use of this system, fringe patterns can be processed and analyzed, and the full field displacement s and strains can be obtained effectively and accura tely. T he system is successfully applied for the residual stress characterization of plain weave textile, shrinkage measurement of concrete material and material property identification of laminate, etc PAGE 16 16 CHAPTER 1 INTRODUCTION Background Acquiring data simultaneously at many points on the specimen is frequently needed for experiments in solid mechanics. A lthough conventional point wise methods such as strain gauges, dial gauges and other mechanical or electrical devices can measure physical parameters s uch as s urface strain and displacement the y can only provide localized measurement. L arge numbers of separate measurements are required to build up an overall picture of the physical parameters; however, when the required number of sensors exceeds 102103 the cost typically becomes prohibitive. B y contrast, moir interferometry, can provide an attractive solution for many applications and can provide fullfield information equivalent to more than 105 independent point wise sensors. W ith moir interferometr y, the data are usually recorded in the form of twodimensional (2D) fringe pattern. S ophisticated digital cameras with high spatial resolution high temporal resolution, and high accuracy have been applied in replace of photographic plates. Moir interferometry is a laser based optical technique that combines the concepts of optical interferometry and geometrical moir I t has been widely applied for studies of composite materials, fracture mechanics, electronic packages, etc. because of its fullfield di splacement measurement capability, high displacement sensitivity, high spatial resolution, and high signal to noise ratio [1] I t is a noncontacting and wholefield method capable of measuring both normal and shear strain. A schematic of the moir interferometry setup can be seen in Figure 11. F igure 12 show s one actual four beam Moir interferometer with a specimen placed in front of it. PAGE 17 17 Figure 11. Four beam Moir Interferometer schematic Figure 12. Four beam Moir Interferometer PAGE 18 18 Fringe patterns recorded via charge coupled device ( CCD) camera are characteristic patterns of light and dark fringes as seen in Figure 13 which represent the intensity of the constructive and destructive interference of light. Mathematically the intensity can be expressed as Equation (11) [2] ,,(,)cos,bmIxyIxyIxyxy ( 11) where (,) xy is the coordinates of the pixel Ixy is the recorde d intensity of each pixel xy ,bIxy is the background intensity ,mIxy is the fringe amplitude, and xy is the phase which is a function of (,) xy xy represent s the fringe order (,) Nxy by ,2(,) xyNxy Figure 13. Typical fringe patterns (scribed circle is 1inch diameter) T he relationship between fringe order and displacement is shown as in Equation (1 2) PAGE 19 19 1 2 1 2x yUxyN f VxyN f ( 12) where Uxy and Vxy are the horizontal and vertical displacement fields. xN and yN are the fringe orders corresponding to the horizontal and vertical displacement fields. f is the specimen grating frequency. U sing the relationships between displacements and engineering strains, the strain can be approximated as shown in Equation (1 3). 1 2 1 2 ,, 11 24x xx y yy y x xyUxy N xy xfx N Vxy xy yfy N UxyVxy N xy yxfyx ( 13) where ,xxxy and ,yyxy are the normal strains in x and y directions, ,xyxy is shear strain x and y are the lengths of gage l ine s (similar to strain gage) in x and y directions. T he traditional technique to analyze a fringe pattern is to manually count the numbers of fringes crossing the gage line and convert it to the displacement and strain information for discrete points [1] A lthough this procedure can be duplicated often to obtain as many data points as possible, it s disadvantages such as low efficiency and lack of fullfield informatio n limit its application. T he term automated strain analysis system refers to the process of converting these fringe patterns into maps representing the parameters of interest: displacement, strain distributions, for example via PAGE 20 20 implementation of mathemat ical algorithms on a digital computer. B ased on the information hidden behind the fringe pattern, there are two basic methods to analyze fringe patterns. O ne is intensity based method and another is phase extraction technique. Literature R eview I n the ear ly days of interferometry, fringe pattern analysis was carried out manually. (Thomas Young [3] was presumably one of the first practitioners of fringe analysis when, in the early 1800s, he measured the spacing of interference fringes in order to calculate the wavelength of light.) During the 1960s, a number of electronic aids [4 5] to fringe analysis were developed. T hese devices allowed a substantial improvement in the accuracy with which fringes could be loc ated within an interferogram but the overall analysis remained essentially manual. D uring the last 30 years, the rapid development of digital image processing equipment has boosted the research on processing fringe pattern automatically. B ased on the principle equation shown in Equation (11), the analysis of a fringe pattern can be achieved through the intensity or phase information. A utomated fringe pattern analysis falls into two categories: one is intensity based and another is based on phase extractio n. M ost early attempts at automated fringe pattern analysis are intensity based [6 28] which is established based on the traditional manual processing method. I t involves detecting and thinning fringe centerlines (fringe skeletons) assigning fringe orders and interpolating fringe orders. A lthough these intensity based techniques are still sometimes used, methods based on phase extraction [29 67] have become more popular. T hey have significant advantages over the intensity based methods [68] : data PAGE 21 21 are obtained over the ful l field, not just at the fringe maxima and minima; the sign of the deformation is given; and immunity from noise is normally better. I ntensity B ased A nalysis T he intensity based fringe analysis method resulted from the traditional manual processing method [1] I n the traditional method, the gage lines are drawn on some points of the fringe pattern and the numbers of fringes crossing these lines are counted separately. Then the displacement and strain in these points can be estimated using Equation (12) and (13). T he matured intensity based procedure includes detecting and thinning fringe skeletons, assigning fringe orders and interpolating fringe orders. A mong them, fringe skeleton detection and fringe order interpolation are the core techniques. T wo kinds of methods can be used to detect fringe skeletons. O ne method involves truncation and binarization of the fringe patterns [21] ; another one is to find the fringe skeleton via detecting the local maxima or minima of fringe intensities [6, 12, 15 20, 28] A lthough many techniques have been proposed and developed for detecting fringe centerlines few studies can be found in the literature for fringe order assignment [10] I t is generally imp ractical to assign fringe orders automatically because the fringe patterns do not contain the fringe order information O ne solution to this problem is to develop a semi automatic method to assign the fringe order which requires the user to detect the ascending or descending directions of the fringe order while applying a small change of load in the experiment [69] Because only the integral fringe orders along the fringe centerlines are assigned, fringe order interpolation is required to obtain fractional fringe orders at every pixel [11] F or fringe order interpolation, onedimension (1D) algorithms based on linear or quadric interpolation are available which are inadequate since the PAGE 22 22 fringe patterns contain twodimensional (2D) information. 2 D fringe order interpolation algorithms are needed to be developed and will be shown in chapter 2. P hase E xtraction M ethods W hen compared to the intensity based methods which only use the information of fringe centerlines phase extraction techniques use the full field fringe information for the analysis. T he phase extraction methods include Fourier transform method [29 41] and phase shifting method [42 67] T he most important advantage of using Fourier transform method to extract the phase information is that only one single fringe pattern is needed for the analysis. H owever, several practical drawbacks limit the application of this method [69] The general Fourier transform method does not work for fringe analysis when the sign of fringe order gradient changes. A lthough ad ding carrier fringe could solve this problem, it increases difficulty in practice because the frequency of the carrier fringe must be controlled accurately. A nother critical limitation of Fourier transform techniques is the lack of capability to handle dis continuities. P hase shifting method is another technique widely used to extract the phase information easily and accurately. F rom the expression of Equation (11), one can see that it is in general impossible to obtain a unique phase distribution from a s ingle fringe pattern because positive phases cannot be distinguished from negative ones without more information. T he near universal solution to this problem is to add to the phase function a known phase ramp which is linear to either time or position. P ha se shifting method uses a series of fringe patterns with phase shifted to calculate the wrapped PAGE 23 23 phase, xy I t is ideal for fringe analysis; however, it also suffers from several problems such as phase unwrapping [70 94] and gradient (strain) calculation. P hase U nwrapping Inconsistent areas from n oise, broken fringes and false fringes exist in most fringe patterns and bring difficulties to phase unwrapping. T o overcome these difficulties, a number of phase unwrapping algorithms have been developed [70 94] in the past 20 years, and new phase unwrapping algorithms are still being proposed for various scientifi c images such as optical shape reconstruction, medical image analysis, geometrical survey, etc. M ost of these algorithms are set for certain problems. Appropriate adoption or improvement of certain algorithms is necessary for t he automated strain analysis system. Itoh s method [70] is the simplest and fastest phase unwrapping procedure. H owever, if any of the pixels are masked or noise exists in the fringe pattern, the process will be interrupted or an error will propagate through subsequent data. O ne class o f algorithm to eliminate the effect from broken fringes is termed branch cut method. These algorithms are based on the identification of the residues [71] in the wrapped phase data and balance them with branch cuts [71 72] T hese algorithms are effective in repairing the broken fringes. A n alternative technique is the so called un weighted least squares phase unwrapping [76, 78] which was developed by the astronomy and synthetic aperture radar signal processing community. T he unwrapped phase map is chosen such that the local phase gradients match the wrapped gradients of the wrapped phase map in a least squ ares sense. T he algorithms above can PAGE 24 24 generally produce a correct solution ; however, they are not able to distinguish false fringes from broken ones since their residues are identical. T he minimum spanning tree method [74, 94] divides the pixels into three groups: 1st group is all the unwrapped pixels (pixels whose phase have been unwrapped) ; 2nd group is the wrapped pixels with at least one neighboring unwrapped pixel; 3rd group is all other wrapped pixels. T he next pixel to be unwrapped is from the 2nd group and should have the minimum phase difference with its unwrapped neighbor among the 1st group. T he preconditionedconjugategradient (PCG) [76, 82] algorithm is the improved unweighted least squares method with the acceleration of the PCG by introducing weights. T he weights are to zeroweight the regions of residues to prevent them from corrupting the unwrapped solution. PCG can offer good performance although it is much slower when compared to branch cut method or minimum spanning tree method. Q uality guided phase unwrapping algorithm [73, 77, 83, 88] does not identify the residues at all. R ather than depending on branch cuts, it relies completely on a quality map [81] to guide the unwrapping paths. The quality map, can determine the consistency of each pixel from the wrapped phase. P hase unwrapping begins at a pixel and grows a region of unwrapped pixels, beginning with the highest quality pixels and ending with the lowest quality pixels. The algorithm is very successful and more efficient than the PCG method, however, it s not guaranteed that the unwrapping paths will not encircle inconsistent areas and thereby introduce spurious discontinuities. So a good quality map is important to this algorithm. And new algorithms to detect and repair these inconsistent areas will be developed to improve the accuracy of the unwrapped PAGE 25 25 phase map. In practice, quality guided phase unwrapping algorithm shows low efficiency for those large fringe patterns. It will be improved in the dissertation to increase its speed of phase unwrapping on large fringe patterns. S ome other algorithms, such as the Flynn s algorithm [80] and the Lpnorm algorithm [79] can also result in good performance for phase unwrapping of fringe patterns. H owever, these methods are less practical and are not discussed here. Strain C alculation I n practice, there is often interest in the strain results which require differentiation of the displacement data (unwrapped phase). G enerally the displacement or unwrapped phase field contains optical and electrical noise. T he noise may not significantly affect the displacement field; however, it can result in large errors in strai n calculation. Finite difference formulae are sensitive to the spacing ( gage length) and normally insufficiently robust. It is necessary to filter the data by, for example, fitting low order polynomials over a small region of the field [95 96] or by Fourier transform low pass filtering [97] T he simplest method is a linear fit over a square subregion [96] A lthough these methods can calculate strain, many discontinuities were observed in practice in the ESA lab Fringe order gradient or strain can also be calculated via differentiating the polynomial function used in the 1D fringe order interpolation [8, 11] 1 D interpolation is effective and fast, however, the result from the interpolations along different directions (x and y direction) are different, which indicate i t is not sufficient to capture the 2 D information. F or simple and uniform fringe patterns, global 2D polynomial interpolation may be used. However, for most real engineering problems, the global 2D polynomial PAGE 26 26 interpolation is not enough to catch the full field phase change due to the complex distribution of the displacement S train calculation can also be done via finite element analysis [98 100] In this method, the displacement fi elds are defined as the boundary condition of the finite element model. H owever, this method can be used only when the material properties are known. I t s accuracy depends on the accuracy of the finite element model. In this dissertation, m ethods based on the global and local surface fit are developed to calculate the strain efficiently with high accuracy Th ey can also be used for fringe order interpolation. O bjectives T he ultimate objective s of the dissertation include the development of the automated st rain analysis system for Moir Interferometry, and the applications of this system for different research projects in the ESA lab. I n order to achieve the goal, the existing algorithms are investigated and improved, and new algorithms are developed. The fo llowing tasks will be completed in this dissertation: Investigation and selection of appropriate algorithms for fringe pattern analysis Improvement of the selected algorithms to increase efficiency and accuracy Experiment al setup for phase shifting and ca libration D evelopment of new algorithms including detection and repair of inconsistent areas in phase map, and strain calculation via global and local surface fit I nvestigation and development of O/DFM method including fringe centerline detection/fringe multiplication/fringe thinning/fringe order assignment/fringe order interpolation/strain calculation from fringe order map PAGE 27 27 D evelopment of an expert software system for fringe pattern analysis Development of a dditional functions to assist fringe pattern analy sis A pplications of this software system in real projects Selected original applications of this system will be completed in this dissertation: Residual stress characterization of plain woven composite Shrinkage measurement of concrete material Material property identification of laminate using full field displacement measurement PAGE 28 28 CHAPTER 2 APPROACH Pre P rocessing of F ringe P attern T he inherent optical noise of an interferogram is shown in Figure 21(A). T he intensities along the red lines in the fring e pattern are also shown in the Figure 21(B) and they are roughly distributed. A lthough the back and white fringe gives good contrast, random noise with large amplitude is evident. T he noise within the fringe pattern will affect the automated fringe analy sis result which is shown in Figure 22. A ppropriate noise filtering algorithms are needed to reduce the effect resulting from the random noise. 0 phase shift 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 A 50 100 150 200 250 300 350 400 450 500 50 100 150 200 250 300 350 400 450 500 50 100 150 200 250 300 350 400 450 500 B Figure 21. Illustration of large random noise existing in fringe pattern. A) Original fringe patt ern B) Intensity distribution along the lines PAGE 29 29 Wrapped phase 3 2 1 0 1 2 3 A 50 100 150 200 250 300 350 400 450 500 2 0 2 50 100 150 200 250 300 350 400 450 500 2 0 2 50 100 150 200 250 300 350 400 450 500 2 0 2 B Figure 22. Illustration of large random noise existing in wrapped phase. A) Wrapped phase. B) Wrapped phase distribution along the lines N oise filtering algorithms fall into two broad categories: spatial domain methods and frequency domain methods [101 102] Spatial domain refers to the image plane itself, and approaches in this category are based on direct manipulation of pixels in an image. F requency domain processing techniques are based on modifying the Fourier transform of an image. I mage P rocessing in the S patial D om ain S patial domain methods are procedures that operate directly on the pixels in the image as shown by the Equation (21). ,, gxyTfxy ( 21) w here fxy is the input image, gxy is the processed image, and T is an operator on fxy defined over some neighborhood of xy Spatial filtering works with the values of the image pixels in the neighborhood and the corresponding values of a subimage (filter, mask, kernel, template, or window) that PAGE 30 30 has the same dimensions as the neighborhood. Figure 23 shows a 3x3 mask 33w and the 3x3 neighborhood of fxy T he process consists simply of moving the filter mask from point to point in an image. A t each point, the response of the filter at that point is calculated using a predefined relationship [101 102] Figure 23. S patial filtering. T he magnified drawing shows a 33 mask and the image section directly under it PAGE 31 31 S moothing linear filters I n general, linear filtering of an image of size MN with a filter mask of size mn is given by the Equation (22): ,, ,ab satb ab satbwstfxsyt gxy wst ( 22) where (1)/2 am and (1)/2 bn wst is the weight coefficient from the mask. T he output of a smoothing linear filter is the weighted average of the pixels within the neighborhood of the filter mask. T hese filters sometimes are called averaging filters. L inear filters can remove the random noise very well. However, it can make the image look blurred. Thus it should be used carefully in practice. O rder statistics filters O rder statistics filters [101102] are nonlinear spatial filters whose response is based on ranking the pixels within the image area encompassed by the filter, and then replacing the value of the center pixel with the value determined by the ranking result. A mong all of the order statistics filters, median filters are quite popular because they provide excellent noise reduction capabilities, with considerably less blurring than linear smoothing filters of similar s ize. M edian filters are particularly effective in the presence of impulse noise, also called salt and pepper noise because of its appearance as white and black dots superimposed on an image. B ut it is not as good as averaging filter for random noise. T he combination of median filter and averaging filter may give the better solution to PAGE 32 32 remove the impulse and random noise. A nd median filters are especially useful when discontinuities exist in the fringe patterns. I mage P rocessing in the F requency D omain L ow frequencies in the Fourier transform are responsible for the general gray level appearance of an image over smooth areas, while high frequencies are responsible for detail, such as edges and noise. T he idea for the noise filtering in frequency domain come s from the removal or attenuation of the frequencies in the frequency domain. A fter the fringe pattern is converted into the frequency domain using Fourier transform as shown in Equation (23), ,, fxyFuv ( 23) T he modification of the frequencies can be done directly using Equation (24), ,,, GuvHuvFuv ( 24) T he multiplication of H and F involves twodimensional functions and is defined on an element by element basis T he filtered image is obtained simply by taking the inverse Fourier transform of Guv I deal low pass filters T he simplest low pass filter is a filter that cuts off all high frequency components of the Fourier transform that are at a distance greater than a specified distance 0D from the origin of the (centered) transform. S uch a filter is called 2D ideal low pass filter (ILPF) and has the transfer function as Equation (25). 0 01 if (,) 0 if DuvD Huv DuvD ( 25) PAGE 33 33 where 0D is a specified nonnegative quantity, and Duv is the distance from point uv to the origin of the frequency rectangle. ILPFs are not very practical because of the severe occurrence of ringing [102] even when only 2% of the total power in the image was removed. T his ringing behavior is a characteristic of ideal filters. Gaussian low pass filters are better to avoid th is phenomena. Gaussian low pass filters (GLPF) T he form of a 2D GLPF is shown in Equation (26). 22 0,/2(,)DuvDHuve ( 26) where 0D is the cutoff frequency, and Duv is the distance from point uv to the origin of the frequency rectangle. W hen 0, DuvD the filter i s down to 0.607 of its maximum value. GLPF can avoid the occurrence of ring well and will be used frequently in practical. S ummary B oth spatial and frequency domain filters are effective in removing the random noise existing in a fringe pattern. B ut they have their own disadvantages. T he wrapped phase obtained from 4 phase shifted fringe patterns which are filtered via 55 averaging filter and GLPF filter separately are shown in Figure 24 and Figure 25. The comparison shows that GLPF filter can always r esult in cleaner wrapped phase map. B ut the frequency domain filter is difficult to apply to those fringe patterns with undefined values. S o the choice of them depends on the specific problem at hand. PAGE 34 34 50 100 150 200 250 300 350 400 450 500 50 100 150 200 250 300 350 400 450 500 3 2 1 0 1 2 3 Figure 24. Wrapped phase map from 4 fringe patterns filtered via 55 averaging filter 50 100 150 200 250 300 350 400 450 500 50 100 150 200 250 300 350 400 450 500 3 2 1 0 1 2 3 Figure 25. Wrapped phase map from 4 fringe patterns filtered via GLPF filter Phase S hifting P hase shifting is a method to extract phase information from fringe patterns. T he following sections will cover aspects of phase shifting. PAGE 35 35 T heory of P hase S hifting T he main elements of a twobeam interferometer are shown schematic ally in Fig ure 2 6. L ight from a coherent light source is divided into the two beams, 1 and 2, which travel along separate paths before being recombined. A photo detector array then measures the resultant intensity distribution. I t is convenient to represent the complex amplitudes of the interfering beams as phasors, as shown in Fig ure 2 7 where 1 a and 2a are the amplitudes and 1 and 2 are the phases of the two beams at a given pixel. Figure 26. Generalized twobeam interferometer Figure 27. Phasor diagram The pixel produces a signal, sxy which is proportional to the intensity of the light falling on to it and can be expressed as Equation (27), 122 22 12 1212, 2cosiisxyaeaeaaaa ( 27) PAGE 36 36 w here is the difference in phase due to the differing optical path lengths encountered by the two beams. F or simplicity will be called phase in the later chapters and the intensity can be mathematically expressed as Equation (28), ,,(,)cos,bmIxyIxyIxyxy ( 28) where (,) xy is the coordinates of the pixel Ixy is the recorded intensity, ,bIxy is the background intensity ,mIxy is the fringe amplitude, and xy is the phase function. F rom the expression of Equation (28), i t s seen that it is general ly impossible to obtain a unique phase distribution from a single fringe pattern. P ositive phase cannot be distinguished from negativ e ones without more information. T he solution to this problem is to add to the phase function a known phase ramp which is linear in either time or position as shown in E quation (2 9) I f the added phase ramp n is linear in time i t i s called temporal phase shifting. I f it i s linear in position, it i s called spatial phase shifting. ,,(,)cos, (0,1,2......1)bmIxyIxyIxyxynnN ( 29) where n is the added phase ramp, N is the total number of fringe patterns with phase shifted, n is the order of the phase shifting pattern. P hase shifting uses a series fringe patterns with shifted phase. A nd the output phase, xy can be calculated using Fourier transform as shown in Equation (210). 1 0 1 0(,)sin(2/) (,) (,)cos(2/)N n n N n nIxynN xyarctg IxynN ( 210) PAGE 37 37 F or example, when 4 N four fringe patterns with shifted phase /2 can be express ed as Equation (211), 0 1 2 3cos, cos, sin, 2 cos,2 cos, 2 cos,3 sin, 2bm bm bm bm bm bm bmIIIxy IIIxyIIxy IIIxyIIxy IIIxyIIxy ( 211) The wrapped phase can then be calculate d as Equation (2 12), 31 02,arctan II xy II ( 212) I n order to use phase shifting, the phase of the fringe pattern needs to be shifted by a certain value. A special stage for our interferometer was designed to add the phase ramp to the fringe pattern. A t the same time, since the phase above is calculated via function arctan, it belongs to and is called wrapped phase. T he wrapped phase cannot be used to calculate displacement. A procedure of phase unwrapping is needed to convert the wrapped phase to the natural unwrapped phase. Experimental S etup The experimental setup for the phase shifting measurement is shown in Fig ure 2 8 It includes a moir interferometer, stage, power supply, specimen grating, CCD camera and computer. The moir interferometer is put on the stage which is controlled. When the stage mov es the fring e pattern is shifted and recorded via CCD camera which is controlled by software in the computer and can record fringe patterns continuously with some time interval. PAGE 38 38 Figure 28. Setup for experiment Stage for P hase S hifting A stage for phase shifting w as designed in the ESA L ab at the University of Florida [103] A s shown in Fig ure 29 the stage is composed of two aluminum plates and four aluminum tubes. T he four tubes were precisely machined and attached to the bottom and top plates at 45o. M agnet ic wi re s of 0.007 diameter with special enamel coating for high temperature were wrapped exactly 200 times around the center of each of the four aluminum tubes then connected together with strain gage wire to complete the circuit. A 0~35 V, 0~5 A adjustable p ower supply, from Pyramid, model PS 32 lab was connected to the circuit. B y turning the power ON and increasing the voltage, the current going through the magnet wire heats the aluminum tubes, which then expands, thus raising the interferometer in a 45o di rection. PAGE 39 39 Figure 29. Stage designed for phase shifting T he stage was calibrated via a time based calibration method. W ith a certain voltage (16V for our research) applied to the stage, the shift ed fringe order and corresponding time was recorded via a high accuracy stop watch program. T he plot of time vs. fringe shift number is illustrated in Figure 210. It showed that a linear shift existed for the low order fringe shifting. T he average time for the first fringe shift is 4.230s. It is divided into se veral intervals depending on how many images w er e taken and used in the phase shifting procedure. 0 1 2 3 4 5 6 7 0 5 10 15 20 25 30 35 40 45 Shift Fringe No.Tims (s)Tims vs. Shift Fringe No. test1 test2 test3 test4 test5 test6 test7 test8 test9 test10 test11 Figure 210. Calibration of the stage PAGE 40 40 Elimination of M iscalibration (N+1 P hase S hifting) E rror from miscalibration of the stage is unavoidable and it mig ht bring significant systematic error to the phase shifting algorithm when the error from miscalibration is large. T his error is called linear phase shift error and several approaches have been developed to reduce the influence of this type of error. A mong these methods, one class of algorithm [104] is especially efficient since it requires only one extra sample (M=N+1) fringe pattern. T his N+1 phase shifting algorithm will be incorporated into the automated strain analysis system. P hase U nwrapping O bviously the wrapped phase obtained from phase shifting belongs to as sho wn in Figure 211. T h e wrapped phase information is not useful in obtaining the displacement information T he procedure to restore wrapped phase back to unwrapped phase is called phase unwrapping. Figure 211. Wrapped phase and unwrapped phase T he relationship between wrapped phase and unwrapped phase can be expressed as Equation (213). A ppropriate multiples of 2 are added to the wrapped phase at each pixel to make the phase map continuously distributed. PAGE 41 41 ,,2,,xyxykxywxy ( 213) where xy is wrapped phase, xy is u nwrapped phase and w is the wrapping operator. P hase unwrapping methods have been extensively studied by numerous researchers The most straight forward is Itoh s method [70] Itoh s method Itoh s method can be easily explained in on e dimension (1D) as follows: w is the wrapping operator (defined by Equation (213)) that wraps the phase into the interval T hus 2, 0,1,...,1, wnnnknnN ( 214) where kn is an array of integers chosen so that n ( 215) T he difference operator is defined as 1, 1, 0,1,...,2, nnn knknknnN ( 216) C omputing the difference of wrapped phases using Equation (214) and (216) yields the Equation (217), 12 wnnnkn ( 217) Again apply the wrapping operator to the above result to obtain, 12 2 wwnwnnknkn ( 218) PAGE 42 42 T he subscripts on 12 and kk are used to distinguish the integ er arrays found by the two wrapping operations. T he result of Equation (218) is the wrapped difference of wrapped phases. B ecause the wrapping operator w produces values in the interval the requirement n ( 219) implies that the term 122 knkn in Equation (2 18) must equal zero. Thus nwwnwn ( 220) which can be manipulated to yield the Equation (221), 1 00m nmwwn ( 221) Equation (221) states that the phase can be unwrapped by integrating the wrapped difference of the wrapped phases. A simple onedimensional phase unwrapping procedure based on Itoh s method is presented as the following algorithm. T his procedure unwraps the phase in the array ,01 iiN Table 21. Procedure of Itohs method of phase unwrapping Step Features Step 1 C ompute the phase differences: 1 Diii for 0,...,2 iN Step 2 C ompute the wrapped phase differences: arctansin,cos i DiDi for 0,...,2 iN Step 3 I nitialize the first unwrapped value: 00 Step 4 U nwrap by summing the wrapped phase differences: 11 iii for 1,...,1 iN PAGE 43 43 Figure 212 illustrates the obtained unwrapped phase using Itoh s method for the pixels along a line in a fringe pattern. O bviously the unwrapped phase is continuous and can be used to calculate displacement and other parameters. 0 20 40 60 80 100 120 140 160 30 20 10 0 10 20 30 40 50 Position (pixel)Phase (Rad)Wrapped phase and unwrapped phase Wrapped phase Unwrapped phase Figure 212. Wrapped phase and unwrapped phase T h e concept above for one dimensional It oh s method can be extended to two dimensional phase unwrapping. T he phase unwrapping algorithm in two dimensions includes a row and a column search of the distribution first by solution of the horizontal discontinuiti es followed by the vertical ones. I toh s method is effective when the fringe patterns are in perfect condition and no broken fringes or false fringes exist. W hen broken or false fringes exit, Itoh s method will result in errors in the unwrapped phase. Fig ure 213 shows the unwrapped phase using Itoh s method for a two dimensional fringe pattern. PAGE 44 44 A B Figure 213. Wrapped phase and unwrapped phase. A) Wrapped phase B) Unwrapped phase obtained using Itoh s method T he inconsistent areas such as broken fringes (D) and false fringes (A,B,C) shown in Figure 212A bring significant errors in the unwrapped phase which is shown in Figure 212B. A better phase unwrapping algorithm is needed to eliminate or limit the affects of these inconsistent areas. A s dis cussed in the literature review chapter, many other phase unwrapping methods, such as Flynn s minimum discontinuity method, cellular automata method [80] minimum spanning tree method [94] and preconditionedconjugategradient (PCG) least square iteration [76, 82] can offer better performance. H owever, these methods are not practical because of their disadvantages which include an unstable result and inefficiency. Quality guided phase unwrapping T he position of inconsistent areas can be determined via detecting and locating the residues. T he calculation of residues is shown in Figure 214 an d Equation (220). PAGE 45 45 Figure 214. Integration direction of residue calculation 1234 1 2 3 4Residue 1,,, 1,11,, ,11,1, ,,1, is the phase wrapping operator. iwhere Wmnmn Wmnmn Wmnmn Wmnmn W ( 222) F or the area without any residues, the integration of the gradient around a closed path should be zero. I f the integration is not zero ( 2 or 2 ), it indicates there is a residue ( negative or positive residue). Figure 215 shows the residues detected using the above method. W hen compared to Figure 213, the position of inconsistent areas and residues accurately agree. PAGE 46 46 Figure 215. Residues detected T he residues will bring wrong 2 or 2 to some areas in the unwrapped phase when using Itoh s method which is shown in Figure 213. An advanced phase unwrapping algorithm is needed to eliminate the inconsistent areas or limit the expansion of the inconsistent areas. T he residues above can tell where the inconsistent areas are. T he information is limited to the inconsistent areas only. A nother parameter, quality map, can determine the consistency of each pixel from the wrapped phase. I n a qualit y map, the areas with low quality represent unreliable phase data and vice versa. T here are many quality maps [81, 93] available such as correlation, pseudocorrelation, phase derivative variance (PDV), etc. Among t hose quality maps, the PDV map can provide the best estimate of consistency of each pixel. T he PDV map is defined by Equation (223), 22 ,,,, xxyy ijmn ijmn mnz kk ( 223) PAGE 47 47 where k is the window size with center pixel mn ,,1,x ijWijij ,1,,y ijWijij x mn and y mn are the averages of these partial derivatives in the kk window. Figure 216B shows the PDV quality map. Obviously in the quality map it is seen that there are several areas with low PDV values c orresponding to the inconsistent areas in the wrapped phase map. W hen phase unwrapping procedure is performed, those pixels with lower PDV values should be unwrapped at the end. Wrapped phase of U field 50 100 150 200 250 300 50 100 150 200 250 300 3 2 1 0 1 2 3 A Quality map of U field (PDV) 50 100 150 200 250 300 50 100 150 200 250 300 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 B Figure 216. Wrapped phase map and its quality map. A) Wrapped phase map B) PDV quality map A fter the PDV quality map is obtained, the quality guided phase unwrapping algorithm can be realized. T he procedure of this algorithm is s ummarized in Table 22 Figure 217 shows a series images which represent the stage of Quali ty guided phase unwrapping procedure. PAGE 48 48 Table 22 Procedure of quality guided phase unwrapping Step Features Step 1 Calculate PDV map mnz and create a binary mask of the same size. S et the value of each point in the mask as 0 (0 means this point is not unwrapped and 1 means it s unwrapped). Step 2 C hoose the pixel with maximum PDV value as the initial pixel and mark the corresponding point in the binary mask as 1, then put information of its neighboring pixels into array ii_ad j, jj_adj (location), and qual_adj (quality value). Step 3 C hoose pixel (ii,jj) from array ii_adj, jj_adj, qual_adj with highest quality within the neighboring pixels, and put its neighboring pixels which has been unwrapped into temporary array ii_adj_tem p, jj_adj_temp, qual_adj_temp. Step 4 D etermine the pixel (ii_base, jj_base) with the highest quality within the temporary array as the base point. Step 5 U nwrap the phase of pixel (ii,jj) based on the neighboring pixel (ii_base, jj_base). Step 6 R epeat steps 3 5 until all the pixels are unwrapped when array qual_adj_temp is empty. 100 Pixels Unwrapped 50 100 150 200 250 300 50 100 150 200 250 300 0 0.5 1 1.5 2 2.5 3 A 1000 Pixels Unwrapped 50 100 150 200 250 300 50 100 150 200 250 300 1 0 1 2 3 B 10000 Pixels Unwrapped 50 100 150 200 250 300 50 100 150 200 250 300 0 5 10 15 20 25 C 50000 Pixels Unwrapped 50 100 150 200 250 300 50 100 150 200 250 300 0 5 10 15 20 25 30 35 40 D 80000 Pixels Unwrapped 50 100 150 200 250 300 50 100 150 200 250 300 10 0 10 20 30 40 E Unwrapped phase 50 100 150 200 250 300 50 100 150 200 250 300 10 0 10 20 30 40 F Figure 217. Process of Quality guided phase unwrapping A) 100 pixels unwrapped. B) 1000 pixels unwrapped. C ) 10000 pixels unwrapped. D ) 50000 pixels unwrapped. E ) 80000 pixels unwrapped. F ) Final unwrapped phase PAGE 49 49 Phase unwrapping starts from the pixel with the highest quality and expands to other pixels according to the order of the quality. W hen compared Figure 217F to Figure 213B, it s seen that the improvement is significant and the affect of inconsistent areas is limited to a small region. I mprove d quality guided phase unwrapping The quality guided phase unwrapping method only relies on the quality map information. T he error from inconsistent areas can be limited to a small area and will not propagate to the whole field. When compared to the Itohs method, the improvement is obvious. Figure 218A and B show the unwrapped phase resulting from the Itohs method and the quality guided phase unwrapping Its seen that the unwrapped phase from the new algorithm is more accurate without streaks emanating from the inconsistent areas with low quality in the quality map. By the use of the new algorithm, the effect of those inconsistent areas is limited and c an be attenuated via interpolation. Unwrapped phase of U field using old PhU 50 100 150 200 250 300 50 100 150 200 250 300 10 0 10 20 30 40 50 60 A Unwrapped phase of U field using Quality Guided PhU 50 100 150 200 250 300 50 100 150 200 250 300 10 0 10 20 30 40 B Figure 21 8 Unw rapped phase map. A) by Itohs method B) by quality guided phase unwrapping PAGE 50 50 The result from the quality guided phase unwrapping method is much better than the Itohs method; however, its ver y time consuming in practic e when applied to a large fringe pattern. T able 2 3 collects the CPU time for the program when applied to fringe pattern with different size. W hen the area is small, the speed of the program is very fast. For the area of 100100 pixels the run time on a PC is less than 5 seconds But for the area of 700700 pixels, the run time increased significantly to about 2.5 hours. The efficiency of the program makes it impractical for large images Table 23 CPU time vs. image size Image size (pixel*pixel) CPU time (s) 100100 4.7344 200200 46.5469 300300 259.4063 400400 843.6094 500500 2149.1 600600 4277.6 700700 8381.3 The solution to this problem is to first divide the large wrapped phase map into small areas, and then apply the quality guided phase unwrapping for each area separately. After all of the areas are unwrapped, they can be connected using boundary pixels. Figure 219A and B show a test on the improved procedure. The 785785 pixel s wrapped phase is divided into small areas and each area is unwrapped separately as shown in Figure 219A Obviously the whole unwrapped phase is not continuous. A procedure of connecting the separate areas is utilized and the final unwrapped phase is shown in Figure 219B PAGE 51 51 Table 24 Improved procedure of quality guided phase unwrapping Step Features Step 0 Divide large wrapped phase map into small areas and do the PhU for each area separately Step 1 Calculate PDV map mnz and create a binary mask of the sam e size. S et the value of each point in the mask as 0 (0 means this point is not unwrapped and 1 means it s unwrapped). Step 2 C hoose the pixel with maximum PDV value as the initial pixel and mark the corresponding point in the binary mask as 1, then put i nformation of its neighboring pixels into array ii_adj, jj_adj (location), and qual_adj (quality value). Step 3 C hoose pixel (ii,jj) from array ii_adj, jj_adj, qual_adj with highest quality within the neighboring pixels, and put its neighboring pixels whi ch has been unwrapped into temporary array ii_adj_temp, jj_adj_temp, qual_adj_temp. Step 4 D etermine the pixel (ii_base, jj_base) with the highest quality within the temporary array as the base point. Step 5 U nwrap the phase of pixel (ii,jj) based on the neighboring pixel (ii_base, jj_base). Step 6 R epeat steps 3 5 until all the pixels are unwrapped when array qual_adj_temp is empty. Step 7 Connect all the separate areas using boundary pixels 100 200 300 400 500 600 700 100 200 300 400 500 600 700 30 20 10 0 10 20 A 100 200 300 400 500 600 700 100 200 300 400 500 600 700 140 120 100 80 60 40 20 0 B Figure 21 9 Unw rapped phase map. A) by Itohs met hod B) by quality guided phase unwrapping The CPU time of unwrapping with for different division sizes is listed in T able 2 5 The improvement of the speed is significant. PAGE 52 52 Table 25 CPU time vs. division size Image size (pixel*pixel) CPU time (s) 7 00 7 00 8381.3 785785 (divided into 150150 areas) 530.2656 785785 (divided into 1 0 0 0 0 areas) 349.3125 Repair of inconsistent areas Some research has been done to repair the broken fringe by placing branch cut s between the positive and negative residues [71 72] However, it is hard to place branch cuts to the right residues since it is hard to detect the right pair of positive and negative residues. And it will also place branch cuts onto the residues from the false fringes. Few studies can be found in literature for false fringe detection and repair. A simple and effective method to detect and repair the false and broken fringes is developed. Figure 220. Removal of false fringe manually A) Wrapped phase wi th false fringe. B) Zoomed in the area with false fringe. C ) Unwrapped phase with big error around the false fringe. D ) Wrapped phase after removal of false fringe. E ) Unwrapped phase with false fringe removed. PAGE 53 53 Quality guided phase unwrapping can restr ict t he error s from inconsistent areas to a small area. Therefore these errors will not propagate to the whole field. However, it cannot attenuate the error from the inconsistent areas such as broken or false fringes. Figure 220A and B show the wrapped phase with a false fringe. With quality guided phase unwrapping, the unwrapped phase of the zoomed area is shown in Figure 220C with phase jump clearly observed. As described previously, the errors caused by the false fringe are limited to the area around. These errors need to be removed or attenuated since they will cause huge errors for gradient (strain) calculation. Eliminat ing the false or broken fringe can be done manually via selection of the area and do the interpolation of the area using the value ar ound the area. Solving L aplaces equation within the area is used to interpolate the phase The mathematical expression Laplaces equation is shown in E quation (2 24) 0, ij is the phase of the points around the area selected. xy is the interpolated phase of the points within area. With the phase of the point in the boundary is given, the phase with distribution of L aplaces equation can be obtained. Figure 220D and E show the wrapped and unwrapped phase with the false fringe area selected and appropriately interpolated. When comparing the unwrapped phase in Figure 220C and E the enhancement resulted from the false fring e elimination procedure is obvious. 22 22 00 ,, xy ijij ( 224) PAGE 54 54 However, this manual selection of the areas is not practical and hard to control the dimension. An automatic method i s developed to detect the inconsistent areas. Detection of inconsistent areas The unwrapped phase for Figure 213 A is shown in Figure 221A with 3 false fringes and 1 broken fringe. The residue map (Figure 213 B) and quality map (Figure 221B) shows the positions of the end points of these false and broken fringes. UnWrapped phase 50 100 150 200 50 100 150 200 20 10 0 10 A PDV quality map 50 100 150 200 50 100 150 200 0.6 0.4 0.2 B Inconsistent area detected 50 100 150 200 50 100 150 200 0 0.5 1 C Dilation of Inconsistent area 50 100 150 200 50 100 150 200 0 0.5 1 D UnWrapped phase after repair 50 100 150 200 50 100 150 200 20 10 0 10 E Difference 50 100 150 200 50 100 150 200 1 0 1 2 F Figure 22 1. Detection and repair of inconsistent areas A) Original unwrapped phase. B) Quality map C ) Inconsistent areas detected. D ) Ext ended inconsistent areas E ) Repaired unwrapped phase map. F ) Difference of unwrapped phase before and after repair PAGE 55 55 The false or broken fringes can be identified by detecting the consistenc y between the neighboring pixels in the unwrapped phase map obtai ned through quality guided phase unwrapping. For the consistent areas, the absolute phase difference between 2 neighboring pixels is very small (<0.3). However, for the inconsistent areas, the absolute phase difference is much larger (>1). This characteris tic can be used to detect the pixels within the inconsistent areas. Starting from the left top corner to the right bottom corner, for each pixel P(i,j), only 3 directions as shown in Figure 222 need to be detected. Figure 222. Detection of inconsist ent areas The inconsistent areas detected are shown in Figure 221C with a threshold value of 1.3. When comparing it with the residue or quality map, the positions of the false or broken fringes agreed well. Additional this method has the ability to detec t those areas with high noise. After the inconsistent areas are detected, phase interpolation using Laplaces equation can be implemented within these areas. However, in order to obtain more trusted boundary pixels for the interpolation, the inconsistent areas can be extended several pixels via image dilation. Figure 221D shows the dilated inconsistent areas with 4 pixels extended. PAGE 56 56 Figure 221E shows the unwrapped phase map after repair of inconsistent areas. The difference of the unwrapped phase before and after repair is shown in Figure 221F. Its seen that the repair s only occur within the inconsistent areas with no affect to other pixels. The p rocedure of detection and repair of inconsistent areas are briefly listed in Table 26 Table 26 Procedure of detection and repair of inconsistent areas Step Features Step 1 Perform quality guided phase unwrapping to obtain the unwrapped phase. And create the mask array with all 0 values. Step 2 Calculate the absolute phase difference between neighboring pixels and compare it with the predefined threshold value. If it is larger than the threshold, these two pixels are marked as inconsistent pixels in the mask array with 1. Step 3 Expand the inconsistent areas in the mask array with dilation. Step 4 Inter polate the phase values within the inconsistent areas using Laplaces equation. The automatic detection and repair of inconsistent areas locates them, remove s the misrepresented phase values and interpolates them with reliable phase information from the boundary pixels. It can reduce the affect caused by the inconsistent areas such as false and broken fringes effectively. Therefore, it can improve the quality of the unwrapped map. In practical application s, the automatic detection method may over detect the inconsistent areas. Currently the solution is to have use interaction when necessary. Future improvement of this method is still needed. PAGE 57 57 S ummary Phase shifting is one of the phase extraction techniques. At least 3 phase shifted fringe patterns are required to calculate the phase information. A stage is designed to add a certain phase ramp to the fringe patterns. PDV quality map can measure the goodness of each pixel. Based on PDV map, quality guided phase unwrapping algorithm is adopted and improved to convert the discontinued wrapped phase to continuous unwrapped phase. The inconsistent areas can be detected. And Laplaces equation is used to interpolate the phase values within the inconsistent areas. Displacement and S train C alculation T he unwrapped phase can be easily converted to displacement via Equation (22 5 ). I t s easily seen that the displacement filed distribution is directly proportional to the unwrapped phase distribution. W hen an in plane displacement field is obtained using Equation (22 5 ), the corresponding strain field is often required to complete the deformation analysis. H owever, the existence of noise always makes it difficult to obtain a perfect strain result. 3 3, 1 10 22 1 10 22U Vxy Uxy m f xy Vxy m f ( 225) T raditional M ethod B ased on the principles of m oir Interferometry the strain distribution can be calculated from the full field dis placement field by Equation (225), PAGE 58 58 6 6 6,, 10 22 ,, 10 22 ,,,, 1 10 2 222U xx V yy UV xyUxyxy xy scale xfx Vxyxy xy scale yfy UxyVxy xyxy xy scale yxfyx ( 226) where is the engineering strain, x and y are the gage lengths (pixels) along x and y direction, scale is the conversion from image coordinates (pixels) to the real specimen dimensio n (mm). G enerally the displacement or unwrapped phase field contains optical and electrical noise. T he noise may not significantly affect the displacement field; however, it can result in large errors when calculating strain. The traditional method to calculate the strain map using gage selection is highly sensitive to the noise existing in the phase map. Gaussian Low Pass Filter may be used to attenuate the noise in the strain map [105] A t the same time, the gage length will also greatly affect the result. A 4 point bending test was conducted here. The material is aluminum and t he dimension of the specimen and the test scheme are shown in Figure 223A. The height of the specimen is 11. 89mm. The fringe pattern taken for analysis is in dimension of 2 2 95 mmx10.40mm. The dimension scale of the fringe pattern i s 41.44pixels/mm. Figure 223B shows 4 consecutive fringe patterns for the U field of the specimen under 4 point bending condition with a line drawn in the position of y=10 pixels PAGE 59 59 A B Figure 223. 4 point bending A) Test s cheme B) 4 consecutive phase shifted fringe patterns The traditional method to calculate strain was used with a gage length of 10 and 30 pixels respectively and strains are illustrated in Figure 224A and B N oise obviously exists in most areas. At the same time, the gage leng th selection affects the result. Generally, large gage length will smooth the result somehow, but still bring great noise into the strain maps as shown in Figure 224C B ased on the result s in Figure 224, a better strain calculation method is needed and t he technique based on surface fit can provide a better solution. PAGE 60 60 Strain xx ( ) 100 200 300 400 500 600 700 800 900 50 100 150 200 250 300 350 400 2500 2000 1500 1000 500 0 500 1000 1500 A Strain xx ( ) 100 200 300 400 500 600 700 800 900 50 100 150 200 250 300 350 400 1500 1000 500 0 500 1000 B 0 100 200 300 400 500 600 700 800 900 200 300 400 500 600 700 800 900 1000 1100 1200 PixelStrain xx ( ) C Figure 224. Calculation of normal strain xx A) Gage length = 1 0 pixels B) Gage length = 3 0 pixels C ) Strain along y=10 pixels PAGE 61 61 Global S urface F it B as ed S train C alculation Generally, the deformation of the specimen measured via m oir Interferometry is continuous, which indicates that the unwrapped phase, displacement, strain and stress should also be continuously distributed. H owever, the existing noise from an optical or electric device will break the continuity and bring significant error in strain calculation. S urface fit can attenuate the noise effectively. A nd most importantly, it can calculate gradient (strain) analytically. T he Equation (22 7 ) sho ws the strain is proportional to the gradient. 6 6 6,, 10 22 ,, 10 22 ,,,, 1 10 2 22U xx V yy UV xyUxy xy xy scale xfx Vxy xy xy scale yfy UxyVxy xyxy xy scale yxfyx ( 227) Some research p apers implement a local surface fit technique to smooth the unwrapped phase field or displacement field. T hese techniques may obtain favorable result over small areas but the overall strain map is not continuous. Global surface fit algorithms are developed here to obtain more accurate strain information. T hin plate splines (TPS) were introduced to geometric design by Duchon [106] and mainly appl ied to the area of computer and vision science [107 109] T he name thin plate spline refers to a physical analogy involving the bending of thin sheet metal. I n the physical setting, the deflection is in the z direct ion, orthogonal to the plane. I n order to apply this idea to the problem of coordinate transformation, one interprets the lifting of the plate as a displacement of the x or y coordinates within the plane. G iven a set of P corresponding points, the TPS warp is described by P +3 parameters which include 3 PAGE 62 62 global affine motion parameters and P coefficients for correspondences of the control points. T hese parameters are computed by solving a linear system, in other words, TPS has closedform sol ution. The theory of thin plate spline is briefly reviewed below, L et iv denote the target function values at locations ii xy in the plane, with 1,2,..., iP T he TPS interpolant xy minimizes the bending energy shown in Eq uation (2 2 8 ) 22222xxxyyy RI dxdy ( 228) A nd has the form 1 1, ,,P xyiii ixyaaxaywUxyxy ( 229) where 2log Urrr I n order for xy to have square integrable second derivatives, it requires that 1 110 0P i i PP iiii iiw wxwy ( 230) T ogether with the interpolation conditions, ,iiixyv this yields a linear system for the TPS coefficients: 0TKPwv POa ( 231) where ,,,ij iiKUxyxy the thi row of P is 1,,iixy O is a 33 matrix of zeros, 0 is a 31 column vector of zeros, w and v are column vectors formed from iw PAGE 63 63 and iv respectively, and a is the column vector with elements 1,,xyaaa D enote the 33 PP matrix by L which is nonsingular [110] and denote the upper left PP block of 1L by 1 pL then it can be shown that 1 TT fpIvLvwKw ( 232) W hen there is n oise in the specified values iv one may wish to relax the exact interpolation requirement by means of regularization. T his is accomplished by minimizing 2 1,n iii iHvxyI ( 233) T he regularization parameter a positive scalar, controls the amount of smoothing; the limiting case of 0 reduces to exact interpolation. A s demonstrated in [111] the TPS coefficients in the regularized case by replacing the matrix K by KI where I is the pp identity matrix. Figure 22 5 shows the normal strain via global surface fit. O bviously it s smoother than that in Figure 224 One limitation of global surface fit is the poor performance on those areas close to the edges perpendicular to the direction of strain calculated. This is due to the poor information extracted in these areas. PAGE 64 64 Strain xx ( ) 100 200 300 400 500 600 700 800 900 50 100 150 200 250 300 350 400 1000 800 600 400 200 0 200 400 600 800 1000 Figure 22 5 N ormal strain xx obtained via global surface fit Global surface fit using TPS is always better to calculate the strain map an d obtain the smooth displacement map. It has a number of advantages: the interpolation is smooth with derivatives of any order; the model has no free parameters that need manual tuning; it has closedform solution for both wrapping and parameter estimation; there is a physical explanation for its energy function. A t the same time, one drawback of the TPS model is that its solution requires the inversion of a large, dense matrix of size P P where P is the number of points in the data set. TPS is an effecti ve tool for the global surface fit and gradient calculation. H owever, the solution requires the inversion of a PP matrix, thus making it impractical for large scale application. O ne obvious solution to this problem is to use the su bset of the data. A trade off between the resolution and computational efficiency has to be made to obtain reasonable result for large images. B spline can also be used to construct the global surface fit. I t is incorporated into the system as one option for the global surface fit. And in order to make use of all the data, local surface fit based strain calculation is developed. PAGE 65 65 Local S urface F it B ased S train C alculation Local surface fit based strain calculation method will first divide the unwrapped phase map into small patches. Surface fit is processed separately on each patch to obtain displacement and strain analytically from the expression of the surface fit. Then all the patches can be connected by the boundary points. Low order polynomials such as bilinear, quadric or cubic function can be used to fit the 2D data within small patches. Here we illustrated the procedure of using quadric function. The format of bilinear and cubic can be found in appendix A and B. The unwrapped phase is assumed to have the quadric distribution within one patch as shown in Equation (234) 22 123456, fijccicjcijcicj ( 234) where fij is the smoothed unwrapped phase. ij is the image coordinate system with unit of pixel. 123456 ,,,,, cccccc are the coefficients calculated from the surface fit. All the points within the patch of size mn pixels are used to form the equation group, 22 111121314115161 22 222122324225262 22 123456, ...... ,mnmnmn mnmnmnmnmnmnfijccicjcijcicj fijccicjcijcicj fijccicjcijcicj ( 235) ,6 mnB 6,1C and ,1 mnY matrix can be formed, 22 111111 22 222222 ,6 221 1 ...... 1mn mnmnmnmnmnmnijijij ijijij B ijijij ( 236) PAGE 66 66 1 2 3 6,1 4 5 6c c c C c c c ( 237) 111 222 ,1, ...... ,mn mnmnmnfij fij Y fij ( 238) Then the coefficients can be calculated as, 1 TT CBBBY ( 239) Finally the gradients (strains) can be calculated analytically as, 346 245 ,2 ,2 ff ijccicj xj ff ijccjci yi ( 240) The negative sign expression of / fy co mes from the opposite direction in image system and rectangular coordinate system. Figure 22 6 shows the normal obtained with the local surface fit with 20x20 patch size and overlap 10. When compared to the result in Figure 22 4 the result from local sur face fit is better and more accurate. When compared to the result in Figure 225, they agree well The result from global surface fit is smoother in most areas. However, local surface fit can obtain more accurate result s on areas close to the edges. This is because local surface can use all the pixels close to the edge to predict the strain. PAGE 67 67 xx 100 200 300 400 500 600 700 800 900 50 100 150 200 250 300 350 400 1000 500 0 500 1000 A xx () 100 200 300 400 500 600 700 800 900 50 100 150 200 250 300 350 400 1000 500 0 500 1000 B Figure 22 6 N ormal strain xx obtained via local surface fit A) Quadric B) Bilinear In real applications of local surface fit met hod, the selection of the function may affect the result a little as shown in Figure 226. All of these functions are incorporated into the system as options. One big problem of using local surface fit is how to connect the patches smoothly. There are big variances of the strain values along the neighboring boundary. A procedure of overlapping between the neighboring patches and then calculating the average values within the overlapped areas is developed to attenuate the variance along the boundary. This procedure works well for most applications. PAGE 68 68 S ummary The traditional method to calculate the strain map is highly sensitive to the noise existing in the phase map. A t the same time, the gage length will also greatly affect the result. The idea of using glob al surface fit technique to smooth the unwrapped phase map can attenuate the noise effectively. A nd most importantly, it can calculate gradient (strain) analytically. TPS was chosen for global surface fit because it is insensitive to noise in the data and it has high capability of constructing complex surface shapes. Theoretically all of the pixels obtained in the CCD camera can be used in TPS surface fit to keep the spatial and strain resolution; however, the calculation of parameters of TPS may be timeconsuming or even fail for large images. O ne obvious solution to this problem is to use the subset of the data, which will reduce spatial resolution somehow. A trade off between the resolution and computational efficiency has to be made to obtain reasonable result s for large images. Local surface fit using polynomial functions is also developed here to utilize all the pixels. Among the polynomial functions, the quadric function has the best performance. A p rocedure of overlapping is developed to reduce the discontinuity of the points along the boundary of each patch. When comparing two methods, the data points along the line of y=10 pixel are extracted and plotted together in Figure 22 7 Since this is a 4point bending test, the strain along y=10 is approxi mate 98 1 ( average strain along the line y=10 calculated manually with 54 fringes within 22.95mm ). The plot shows both methods can obtain more reliable strain values than traditional method as shown in Figure 224C for most areas w ith strain values close to the theoretical strain. The ranges of strain for most areas from global and local surface fit are 98130 and 98160 respectively, which PAGE 69 69 indicates that global surface fit can give smoother or more accurate strain calculation than the local surface fit method for most areas. However, local surface fit has better performance on those areas close to the edges. At the same time, the results indicate moir has a better resolution of strain than other optical method such as digital image correlation which generally has strain resolution of 200 0 100 200 300 400 500 600 700 800 900 200 300 400 500 600 700 800 900 1000 1100 1200 X: 657 Y: 921.1 Pixel X: 394 Y: 1043 X: 421 Y: 1017 X: 627 Y: 961.3 Strain xx ( ) Theoretical Strain Global Surface Fit Local Surface Fit Figure 22 7 Plot of strain along y=10 pixels Optical/D igital F ringe M ultiplication (O/DFM) M ethod T he following sections will c over the procedure associated with the O/DFM method to analyze the fringe patterns. T heory of O/DFM M ethod O/DFM is an improved intensity based technique to analyze fringe patterns. I n O/DFM, a series of N fringe patterns with phase shift ing are used [1, 21] F or a series of N shifted fringe patterns, the phase of each m oir pattern is shifted by 2/ N relative to its neighbors. W hen N is an even number, the fringe patterns can be divided into two groups: the patterns of the first half and their complements. The intensity distributions of these two groups are shown in Equation (241 ), PAGE 70 70 2 ,,(,)cos, (0,1,2......1) 2 2 ,,(,)cos,ibm ibmi IxyIxyIxyxy N N i i IxyIxyIxyxy N ( 241) where iI is the intensity distribution of the thi shifted pattern which is shifted by 2/ iN with respect to the original pattern. iI is the intensity distribution of the corresponding complementary pattern which is shifted by wit h respect to the thi shifted pattern. 0 2 4 6 8 10 12 14 16 3 4 5 2 Moire fringes with phase shifted y1=4+cos(t),y2=4+cos(t+) 0 2 4 6 8 10 12 14 16 2 0 2 Substract shited intensity y3=y1y2=2cos(t) 0 2 4 6 8 10 12 14 16 0 1 2 Take absolute values y4=abs(y3)=2cos(t) 0 2 4 6 8 10 12 14 16 0 0.5 1 Binarize near y5=0 Figure 22 8 Explanation of O/DFM A) 2 complementary fringes B) Subtraction of each other C ) Absolute value of B) D ) After binarizing A B C D PAGE 71 71 T he complementary patterns are used in O/DF M method and explained in Figure 2 2 8 W hen the two complementary patterns are subtracted from each other at each point, as illustrated in Figure 22 8 B. Then it takes the absolute values as illustrated in Figure 22 8 C. The algorithm proceeds by truncating the data near 0 and binarizing intensities of 0 and 1 to points below and above the truncation value [13 14] respectively, as shown in Figure 22 8 D. Obviously the result is a sharpened contour map that has twice as many contour lines as the number of fringes in the initial pattern. T he sharpened contours occur at the crossing points of the complementary fringe patterns, where the intensities are equal. O/DFM method can provide N times fring e density than the original fringe pattern which will bring more accurate analysis result. T he O/DFM algorithm in this dissertation includes the following 5 steps; fringe skeleton detection, fringe thinning, fringe order assignment, fringe order interpolat ion, and strain calculation. 4 consecutive fringe patterns with /2 phase ramp are used below for illustration. F ringe S keleton D etection B ecause of the existence of noise in the fringe pattern, the binarizing procedure shown in Fi gure 22 8 generally does not work well. Figure 22 9 shows the skeleton detected using binarizing method for only 2 complementary fringe patterns 01II O bviously, the skeleton thickness is not uniform and in many places it is very thic k. L ocal minimum detection algorithm [6, 12, 1520, 28] can provide a better solution to the fringe skeleton detection. I n this algorithm, the fringe skeleton is represented as the local minimum pixels which can be found using the following procedure. PAGE 72 72 Fringe skeleton after truncated with 0.1 100 200 300 400 500 600 700 800 100 200 300 400 500 600 700 800 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Figure 22 9 Fringe skeleton detected using binarizing method for 02II Figure 230. Fringe skeleton detected local minimum PAGE 73 73 F igure 230 shows the four directions which will be tested. ,2,1,0,1,2ijPij represents the intensity of the pixels around the center pixel 0,0P. T he test of each direction is shown in Equation (240 ) I f 2 of the above 4 conditions are satisfied, pixel 0,0 P is r egarded as the local minimum pixel and its value is set to 0. 1,00,01,01,20,21,2 1,00,01,01,20,21,2 0,10,00,12,12,02,1 0,10,00,12,12,02,1 1,10,01,12, direction: direction: direction: PPPPPP xPPPPPP PPPPPP y PPPPPP PPPP xy 12,21,2 1,10,01,12,12,21,2 1,10,01,12,12,21,2 1,10,01,12,12,21,2 direction: PP PPPPPP PPPPPP xy PPPPPP ( 242) Fringe skelelton of U field for 2 fringes 100 200 300 400 500 600 700 800 100 200 300 400 500 600 700 800 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Figure 231. Fringe skeleton detected using local minimum detection method for 02II PAGE 74 74 Figure 231 shows the fringe skeleton resulted from this method and obviously the skeleton is much thinner than that from binarizing procedure as shown in Figure 22 9 F igure 231 also shows some isolate noise and the thickness of the skeleton is more than 1 pixel in most fringes. A ppropriate fringe thinning algorithm is needed to obtain the best fringe skeleton. F ringe T hinning G enerally, the fringe skeleton obtained from the binarizing or local minimum detection method is wider than the width of one pixel. Fringe thinning is the process to reduce the skeleton width and thus to obtain a true fringe skeleton [22 23, 112113] Figure 232 shows the fringe skeleton after thinning. W hen compa red to Fig. 2 31 it s thinner and cleaner. Figure 232. Fringe skeleton after fringe thinning for 02II PAGE 75 75 F ringe O rder A ssignment F ringe skeletons represent the contours of equal displacements and they don t contain any informatio n of fringe orders. T herefore, fringe order assignment is needed for the wholefield fractional fringe order calculation. S emi automatic fringe order assignment algorithms were developed by some researchers [10 11] T hese methods require the user to inform the program the tilting direction of the fringe skeleton interactively. E xtra ef fort is needed to determine the tilting direction beside the phase shifting experiment. A new method to assign the fringe orders is developed here. I t use s the unwrapped phase obtained from phase shifting. S ince the unwrapped phase map cont ains the full f ield information it also contains the contours of equal displacement. 100 200 300 400 500 600 700 800 100 200 300 400 500 600 700 800 10 8 6 4 2 0 2 4 6 8 10 Figure 233. Fringe order assignment PAGE 76 76 F ringe O rder I nterpolation F ringe order assignment procedure can only determine the fringe order of the skeletons. T o obtain the full field fra ctional fringe orders, appropriate fringe order interpolation is needed. A 1 D interpolation algorithm using cubic spline has been applied to fringe order interpolations [8, 11] I t can result in continuous 1st and 2nd derivative of the fringe order. H owever, it s not s ufficient to describe the 2D fringe patterns. T he main disadvantage of 1 D interpolation is that interpolations along x direction and y direction do not yield the exactly same results. S o a 2D fringe order interpolation is necessary. Displacement: U field (m) 100 200 300 400 500 600 700 800 100 200 300 400 500 600 700 800 4 3 2 1 0 1 2 3 4 Figure 234. Fr inge order interpolation T he idea of global or local surface fit from strain calculation of phase shifting is reused here. T he surface fit interpolation can obtain the full field information, and provide an effective way to calculate the gradient (strain) of the fringe order. PAGE 77 77 S train C alculation B ased on the fringe interpolation above using surface fit method, the strain (gradient) can be easily obtained as shown in Figure 235B Strain xx ( ) 100 200 300 400 500 600 700 800 100 200 300 400 500 600 700 800 440 420 400 380 360 340 320 A Strain xx ( ) 100 200 300 400 500 600 700 800 100 200 300 400 500 600 700 800 440 420 400 380 360 340 320 B Figure 235. Strain calculated. A) By phase shifting. B) By O/DFM PAGE 78 78 Th e strain map is smoother than that obtained by phase shifting as shown in Figure 235A. This is because most data points in O/DFM method are from the interpolation. Fewer seed points can be used to construct global surface fit. S ummary T he theory of O/DFM to obtain high density fringe skeleton is reviewed. A ll of the necessary procedures associated with this technique are developed. I t requires even number of fringe patterns with phase shifted. S ince it uses fringe order interpolation to calculate displacement and strain, the effect of noise can attenuate H owever, to determine the fringe order accurately and automatically, it requires the unwrapped phase map obtained from phase shifting method. This way the O/DFM method is a hybrid method combining the tec hniques of phase extraction and intensity based methods And since most global surface interpolation is used to obtain the fringe order for most of the data points, smoother displacement and strain field can be obtained. S ummary This chapter investigated the existing fringe analysis techniques and improved or developed new techniques to analyze the fringe patterns automatically. T he following subjects are coved in this chapter, Pre processing of fringe pattern P hase shifting Displacement and strain calcu lation O/DFM Pre processing of fringe pattern involves spatial or frequency filters applied on fringe patterns to remove various noise. P hase shifting is a one phase extraction PAGE 79 79 method to analyze fringe pattern automatically. I t requires at least 3 consecut ive fringe patterns with phase shifted. A stage was designed to apply phase ramp to the fringe patterns. A phase unwrapping algorithm based on quality map was adopted and improved to calculate the unwrapped phase map. Inconsistent areas in the unwrapped phase map were detected and repaired. G lobal surface fit and local surface fit based on thin plate theory and polynomials separately w ere developed to calcul ate displacement and strain maps with higher accuracy instead of the qualitative analysis provided by most of the existing systems O/DFM method is a hybrid method suitable to analyze fringe patterns with low density of fringes. B y combining with the unwrapped phase result from phase shifting, O/DFM can also analyze fringe pattern automatically. PAGE 80 80 CHAPTER 3 AUTOMATED STRAIN ANA LYSIS SYSTEM P latform T he automated strain analysis system is developed as a form of software. T he software is developed using MATLAB R200 9 a. It is based on Windows Graphic User Interface (GUI) which is menudriven and mouse driven, and can be executed on Microsoft Windows NT/2K/XP /Win 7 operating systems. F unctions of the S ystem T he software system includes all the algorithms discussed in previous chapters and other additional image processing functions. T he main functions of the sof tware are listed in Table 31, Table 31. List of f unctions in the software Category Features File U field (N) V field (N) U field (N+1) V field (N+1) U field (O/DFM) V field (O/DFM) Close Image Processing Noise Smooth Image Crop Boundary No Bound ary Define Boundary of Circle/Ellipse/Rectangle Define Boundary using Line Define Boundary using Spline Load Boundary Phase Unwrapping Quality Guided PU Phase Smoothing O/DFM Fringe Skeleton Fringe Thinning Fringe Order Assignment Fringe Orde r Interpolation Tool Data Along A Line Phase Repair Binary Image Repair Complex Boundary PAGE 81 81 Displacement/Strain U Displacement V Displacement Strain xx Strain yy dU/dy dV/dx Shear Strain R Displacement Theta Displacement Strain rr Strain Theta S hear Strain rTheta Display Filtered Fringe Wrapped Phase Unwrapped Phase U Displacement V Displacement Deformation Plot Strain xx Strain yy Shear Strain exy R Displacement Theta Displacement Strain rr Strain Theta Shear Strain rTheta Help About Demo vi deo D emonstration of the S ystem A demonstration of the system for the analysis of the fringe patterns from a tension test on an openhole aluminum specimen is shown here. The dimensions of the test specimen and fringe patterns of the V field are shown in Figure 31. PAGE 82 82 Figure 31. Open hole test specimen and fringe patterns Figure 32 ~8 show some screen captures of the software. Figure 32 Open image files PAGE 83 83 Figure 33 Coordinate system translating parameters Figure 34 Noise removal PAGE 84 84 Fi gure 35 Define the boundary PAGE 85 85 A B Figure 36 Quality guided phase unwrapping A) Unwrapping progress B) Final unwrapped phase PAGE 86 86 Figure 37 Phase smoothing PAGE 87 87 A 0 50 100 150 200 250 300 350 400 1000 1500 2000 2500 3000 3500 4000 4500 X: 364 Y: 4021 PixelStrain xx ( ) X: 1 Y: 1446 X: 114 Y: 1558 B Figure 38. Data extraction. A) Define the data path. B) Data along the data path PAGE 88 88 Validation Two tests, 4 point bending test and openhole test, were conducted to validate the system. The fringe patterns for the tests were analyzed using the system. The results are compared with those results obtained from manual calculation or finite element (FE) simulation. 4 p oint Bending Test The dimension of the 4 point bending test specimen is shown in Figure 2 23. The fringe patterns for the U and V field are shown in Figure 39. A B Figure 39 Fringe pattern. A) U field. B) V fie ld These fringe patterns were inputted into the automated strain analysis system and the result strain maps are shown in Figure 310. PAGE 89 89 Strain xx ( ) X: 475 Y: 1 Index: 1062 RGB: 0, 0, 0.547 X: 475 Y: 431 Index: 1028 RGB: 0.547, 0, 0 100 200 300 400 500 600 700 800 900 50 100 150 200 250 300 350 400 1000 500 0 500 1000 A X: 475 Y: 1 Index: 111.4 RGB: 0.828, 0, 0 Strain yy ( ) X: 476 Y: 431 Index: 830.9 RGB: 0, 0, 0.703 X: 1 Y: 216 Index: 40.25 RGB: 1, 0.0938, 0 X: 951 Y: 216 Index: 95.91 RGB: 0.891, 0, 0 100 200 300 400 500 600 700 800 900 50 100 150 200 250 300 350 400 800 600 400 200 0 200 B X: 169 Y: 45 Index: 105.4 RGB: 1, 0.75, 0 Shear Strain xy ( ) X: 924 Y: 48 Index: 4.995 RGB: 0.281, 1, 0.719 X: 577 Y: 33 Index: 41.4 RGB: 0.688, 1, 0.313 X: 445 Y: 408 Index: 30.46 RGB: 0.594, 1, 0.406 X: 71 Y: 397 Index: 47.86 RGB: 0, 0.891, 1 X: 917 Y: 361 Index: 5.898 RGB: 0.375, 1, 0.625 X: 460 Y: 258 Index: 41.42 RGB: 0.688, 1, 0.313 X: 22 Y: 258 Index: 105.8 RGB: 1, 0.734, 0 100 200 300 400 500 600 700 800 900 50 100 150 200 250 300 350 400 200 100 0 100 200 C Figure 310. Strain maps from automated strain analysis system. A) Strain xx B) St rain yy C) Strain xy PAGE 90 90 A FE model was developed to simulate the 4point bending test. The strain maps within the same area are shown in Figure 311. When compared to the results in Figure 3 10, they agree very well for most areas. A B C Figure 311. Strain maps from FE A) Strain xx B) Strain yy C) Strain xy Furthermore, the result can be compared to the manual calculation for se veral discrete points as shown in Figure 312. PAGE 91 91 A X: 200 Y: 75 Index: 675.1 RGB: 0, 0.25, 1 Strain xx ( ) X: 200 Y: 325 Index: 535.5 RGB: 1, 0.484, 0 X: 478 Y: 50 Index: 834 RGB: 0, 0, 0.953 X: 481 Y: 400 Index: 891.1 RGB: 0.813, 0, 0 X: 484 Y: 2 Index: 1061 RGB: 0, 0, 0.516 X: 480 Y: 430 Index: 1021 RGB: 0.563, 0, 0 X: 800 Y: 70 Index: 687.9 RGB: 0, 0.219, 1 X: 800 Y: 385 Index: 819.6 RGB: 0.938, 0, 0 100 200 300 400 500 600 700 800 900 50 100 150 200 250 300 350 400 1000 500 0 500 1000 B Figure 312. Strain result from manual calculation and analysis system A) Strain xx from manual calculation. B) Strainxx from the analysis system Th e strain results from manual calculation and the automated strain analysis system match very well for these discrete points. Openhole Tension Test The dimension of the 4 point bending test specimen is shown in Figure 31. The fringe pattern for the V field is also shown in Figure 31. These fringe patterns were inputted into the automated strain analysis system and the result strain map ( yy ) are shown in Figure 31 3A A FE model was developed to PAGE 92 92 simulate the Openhole tension tes t. The strain map ( yy ) within the same area are shown in Figure 313B. It matches very well with the result from the analysis system for most areas except the extremely high strains in the small areas close to the edge of the hole due to the low resolution and low quality of fringe patterns within these areas 100 200 300 400 500 600 700 800 900 50 100 150 200 250 300 350 400 450 500 0 500 1000 1500 2000 2500 3000 3500 4000 A B Figure 31 3 Strain map ( yy ) from automated strain analysis system and FE simulation A) From automated strain analysis system B) From FE simulation PAGE 93 93 Furthermore, the result can be compared to the manual calculation for several discrete points as shown in Figure 31 4 A X: 1 Y: 256 Index: 1446 RGB: 0.0938, 1, 0.906 X: 51 Y: 256 Index: 1511 RGB: 0.156, 1, 0.844 X: 101 Y: 256 Index: 1550 RGB: 0.188, 1, 0.813 X: 151 Y: 256 Index: 1590 RGB: 0.219, 1, 0.781 X: 201 Y: 256 Index: 1621 RGB: 0.25, 1, 0.75 X: 251 Y: 256 Index: 1705 RGB: 0.328, 1, 0.672 X: 301 Y: 256 Index: 1902 RGB: 0.516, 1, 0.484 X: 331 Y: 256 Index: 2543 RGB: 1, 0.906, 0 X: 351 Y: 256 Index: 3249 RGB: 1, 0.25, 0 X: 364 Y: 256 Index: 4021 RGB: 0.547, 0, 0 X: 454 Y: 50 Index: 683.1 RGB: 0, 0.391, 1 X: 200 Y: 100 Index: 1710 RGB: 0.328, 1, 0.672 B Figure 31 4 Strain result from manual calculation and analysis system. A) Strain xx from manual calculation. B) Strainxx from the analysis system PAGE 94 94 The strain results from manual calculation and the automated strain analysis system match very well ( 5% variance) for these discrete points The strains along the narrowest section from the automated strain analysis system, manual calculation and FE simulation are shown in Figure 315. Good agreement between them is observed. 0 1 2 3 4 5 6 7 8 9 10 11 12 13 1000 1500 2000 2500 3000 3500 4000 4500 X(mm)Strain yy ( ) Phase shifting Strain gage FE simulation Figure 315. Strains along the narrowest section T he stress or strain concentration factors fro m the FE simulation, manual calculation and the analysis system are calculated and shown in Table 32 They agree well with each other. Table 32 Concentration for the openhole aluminum test Theoretical (Mechanics of Materials) FE Manual System Stres s Concentration 2.65 2. 54 Strain Concentration 2. 64 2. 6 2 2. 59 S ummary The automated strain analysis system is developed as a Windows GUI based software to analyze fringe patterns automatically It includes all the algorithms PAGE 95 95 presented in the previous chapters. With the use of this expert system, full field displacement and strain information can be effectively and accurately extracted from the phase shifted fringe patterns. And the system can be updated step by step with the development or improvement of new algorithms in the future. Two tests were conducted to validate the system. The fringe patterns were analyzed via the system and manual calculation. FE models were conducted for both tests. The results obtained from the system, manual calculation a nd FE simulation are compared and they show great agreement with each other. PAGE 96 96 CHAPTER 4 APPLICATIONS OF AUTO ATED STRAIN ANALYSIS SYSTEM R esidual S tress C haracterization of P lain W oven C omposites Introduction Composite materials are widely used by various industries. T here are two categories of composites, plies composed of unidirectional or randomly scattered chopped fibers (laminates), and fiber bundles woven into a regular pattern (textiles). A s the desire increases for advanced material properties tha t composites provide, more research needs to be performed so that they can be fully characterized and understood. T he residual stresses resulting from the manufacturing process of composite gain a lot of attention because they will reduce the life time of the composite structures and neglecting them will overestimate the performance prediction. Much research has been devoted to determining the residual strains that develop after the cure of laminate d composites. However, with the growing popularity of woven composites, attention in research must be shifted. In order to fully understand how the woven structure will behave after cure, the residual strains associated with the curing process must be measured. By understanding the residual stresses that develop, one can better assess the structural life remaining in the specimen in question. R esidual stresses within composites come from the requirement of high temperature heating cycle during the manufacturing process. T hose high temperatures are required so t hat the resin can fully heat and complete its polymerization process as well as wet the fibers and finally cure into a hard structure. This necessary process introduces a significant amount of residual stresses into the composite specimen after the heating cycle completes and the specimen cools down to room temperature. T he PAGE 97 97 residual stresses in composites arise from both chemical and thermal shrinkage [114] T he chemical shrinkage is due to polymerization of the resin in which the two monomers in the epoxy come together to form the final compound. T he majority of the chemical shrinkage occurs during the initial heat ing period of the curing cycle and therefore does not contribute to the overall residual stresses. T his is because stresses cannot exist before the composite fully cures and the resin transforms to the solid state. H owever, there is some additional chemica l shrinkage that does contribute to the final residual stress levels in the specimen that occurs after the solidification of the resin. O n the ply level, a significant amount of the residual stresses develop because of the mismatch that naturally exists between the coefficients of thermal expansion (CTE) of the two materials. T he resin typically has a much higher CTE value than the fiber; therefore a large amount of thermally induced shrinkage would develop in the matrix. H owever, the higher stiffness fiber s restrict the contraction of the matrix. T hat restriction is what introduces the residual stresses into the composite. U nlike the uniform strain distribution of laminates, significant strain variation within the representative volume element (RVE) of tex tile composites, as shown in Figure 41, is proven under a tension loading test [1, 105] A nd the residual strain distribution should have the same difference between laminate and woven textile. S ince the strain var iation exists in woven textiles, an experimental technique with full field measurement capability is desired. PAGE 98 98 Figure 41. Representative volume element of the plain woven textile M any methods exist to determine the residual str esses of composite mater ials. T hey can be divided into destructive methods and nondestructive methods. D estructive methods include holedrilling [115116] and first ply failure test [117 118] N ondestructive methods include X Ray diffraction [119] imbedding sensors [120] neutral diffraction method. T he destructive methods will destroy the specimen and can only do the measurement of several points. X Ray diffraction wil l introduce foreign materials to the specimen. T he result from imbedding sensors is localized and the expense of neutron diffraction method is extremely high. I n order to measure the residual strain for woven textiles, a nondestructive method with full fi eld measurement capability is desired. C ure reference method (CRM) is such a method. Classical laminate theory (CLT) is suitable to calculate the residual stresses of traditional laminate composites once the material properties and strain data are known [103, 121] However, it cannot be applied to woven composite because the complexity of woven geometry no longer allows these straightforward calculations To predict the residual stress es associated with the residual strains of woven composite, Finite PAGE 99 99 Element (FE) mod eling is the preferred analysis method [122 123] And FE model can compare the residual strains obtained experimentally by simulating the curing cycle. Cure R eference M ethod CRM was originally developed by Ifju, et al. at the University of Florida as a nondestructive technique to measure full field residual strains of composite materials as the composite cools via m oir i nterferometry [124] By using CRM, the residual strains could be measured for the woven composite geometry on a full field scale, and for the RVE Because of the regular repeating pattern exhibited by woven composites, it has been shown that a distinct regular fringe pattern should develop that would aid in the strain characterization [125] Additionally, the use of CRM does not reinforce the specimen and avoids destruction of the specimen. The focus of this work is on plain woven composites. The prepreg is CYCOM 97714A made by CYTEC company. T o use moir interferometry in CRM, the autoclave tool shown in Figure 42 was used firstly to tune the interferometer to obtain null field fringe pattern. T hen it was removed and the specimen was placed in the same position. B y this way, the absolute residual strain in room temperature referred to the cure temperature can be measured. Figure 42. C ure reference method PAGE 100 100 T he displacement field is obtained via CRM which requires a moir interferometry diffraction grating be attached to the composite during cure. A diffraction grating was transferred from the autoclave tool to the prepreg at the cure temper ature while the specimen was at a free stress state. When the specimen cools to room temperature, the grating records the complete strain development. Ifju et al. [124] described a procedure of creat ing an autoclave tool for the CRM Each grating made by the process required approximately 78 hours includ ing the cure times, alumini zation steps and dwelling times. S ignificant w ork was done by Strickl and [126] to reduce the production time to approximately 45 hours. This reduction allow ed quicker production and enhanced the grating quality with fewer processes and fewer opportunities for errors as shown in Figure 4 3 T he technique yielded mixed results in terms of grating quality. Figure 43 4 CRM gratings PAGE 101 101 Result of R esidual S trains from E xperiment O ne of the specimens in Figure 43 was placed in front of the moire interferometer and fringe patterns are reco rded by the CCD camera and shown in Fig. 44. A B Figure 44. Fringe patterns. A) U field B) V field The automated strain analysis system was used to obtain the displacement and strain maps. Figure 4 5 shows the displacement distribution for the U and V field. U field (m ) 100 200 300 400 500 600 700 800 900 1000 100 200 300 400 500 600 700 800 900 000 20 15 10 5 0 A V field (m) 100 200 300 400 500 600 700 800 900 1000 100 200 300 400 500 600 700 800 900 1000 0 5 10 15 20 B Figure 45 Displacement field. A) U field B) V field PAGE 102 102 Based on the displacement field, the strain field can be obtained as shown in Figure 46 I t s seen that the maximum repeating tensile strains developed were approximately 5 0 0 and the maximum repeating compressive strains were approximately 1400 in both the U and V fields. T he strain map that developed can be explained by relating the trends to the geometry of the woven textile. T here are two main regions within the RVE, the resin and the fiber regions. T he resin zone expands the entire length of the RVE between fiber bundles. H owever, along that length, the resin depth from the surface changes as it crosses over transverse bundles. T he fiber regions have much lower surface resin content levels since they are much closer to the surface. Strain xx ( ) 100 200 300 400 500 600 100 200 300 400 500 600 1400 1200 1000 800 600 400 200 0 200 400 A Strain yy ( ) 100 200 300 400 500 600 100 200 300 400 500 600 1400 1200 1000 800 600 400 200 0 200 400 B Figure 46 Strain map. A) Strain xx B) Strain yy T he regions of tensile strains oc curred along the top of the fiber bundles in the fiber regions in both fields of view, whereas the compressive strains were developing throughout the resin zones. T he CTE value for the resin is significantly larger than that of the fiber, so as the composi te cools to room temperature, a considerable amount of PAGE 103 103 residual strains develop within the resin because the fiber bundles restrict free contraction in those areas. T he tensile regions that developed along the fiber bundles were a direct result from the co mpressive resin regions. T o maintain static equilibrium, the thin resin zones that cover the fibers must undergo tension to accommodate for the large compressive areas. A t the same time, it can also be noted that the strain patterns and values are similar for the U and V fields. T hat was expected to occur because of the symmetry of the woven geometry. The strain within one RVE can be extracted via averaging the four RVEs in Figure 4 6 and they are shown in Figure 47 Furthermore, the average strain of the RVE can be obtained as 596 for U field and 450 for V field. Strain xx ( ) 50 100 150 200 250 300 50 100 150 200 250 300 1200 1000 800 600 400 200 0 200 A Strain yy ( ) 50 100 150 200 250 300 50 100 150 200 250 300 1000 800 600 400 200 0 B Figure 47 Strain map within one RVE A) Strain xx B) Strain yy Finite E lement M odel D escription T o predict the corresponding residual stresses, finite element model can be used to simulate the complex geometry of textile composites. The finite element model developed by Karkkainen [123] of the representative volume element for the plain weave composite was used with the adjustment o f the dimension and boundary PAGE 104 104 conditions. The geometry shape of the FE model for the representative volume element (RVE) as shown in Fig ure 4 8 was taken from a literature source that had documented the dimensions required when constructing the woven patter n [127] The real dimensions were measured via digital caliper along with the assistance of microscope. 0.1cos2 0,4() 4 1(), 0.05() 0.1(), 0.36()ab fmx z xmm Rmm Rmm hmmhmm Figure 4 8 Geometry shape Parametric modeling using Pro/Engineering is used to construct the 3D structure, which can better simulate the geometry of textile and easily adjust the geometry which are shown in Figure 49 By this way, the dimension of the geomet ry structure can be easily modified to make it close to the real dimension. T he CAD file is imported into Abaqus software separately and assembled together. Contact condition between fiber and matrix is defined as glued together under the assumption that n o slip occur s between them. PAGE 105 105 A B Figure 4 9 Geometry model of fiber and matrix A) Fiber B) Matrix Two materials were created within the model for the matrix and the yarns The yarn was assigned the properties of a unidirectional laminate, AS4/35016 taken from experimental data [103] However, because of the undulation that occurs in the weave, the properties could not be assigned directly to the entire fiber bundle. A local material coordinate system to each element was built to guarantee the laminat e properties aligned correctly with the local direction of the fibers. Only one material needed to be defined for the weave by this mean. The temperature dependent properties for the matrix, 35016, were taken from the literature [128] To simulate the curing cycle, the application of temperature and pressure fields w ere required to simulate both the change in temperature and releasing of the vacuum surrounding the specimen. When the temperature dependent Coefficient of Thermal Expansion (CTE) data was input into the model, the reference temperature was set to be 1 26.7oC so that any temperature field applied would simulate cooling from the cure temperature. The room temperature conditions were modeled by applying a full field temperature of 22oC and a hydrostatic pressure load of 101.3kPa. To avoid rigid body motion, while allowing full deformations of the model, the center node was fixed in all PAGE 106 106 three translational directions. All of the material properties used in the model can be found in Table 4 1. Table 41. Material properties used in FE model for plain weave composite RVE Temperature ( C ) 1 () EGPa 2 () EGPa 12() GGPa 12 1(/) C 2(/) C AS4/3 501 6 120 138.00 7.67 4.07 0.30 0.59 2 7 32 20 138.00 9.3 4 5.30 0.30 0.45 23.03 35016 120 3. 1 0 0. 41 52.94 52.94 20 4 .5 5 0.3 6 41.30 41.30 The automated strain analysis system provides displacement and strain information at every pixel. Theoretically the full field displacement data can be input to the FE model to obtain the resulting full field stress data. How ever, due to the large number of data points from the experiment and number of nodes in the FE model, the Periodic Boundary Conditions (PBC) was used to import the displacement information as follows: 1 2, ,RL RL TB TBuu vv uu vv ( 41) where R,L,T,B and 1 and 2 are the right, left, top, bottom edges and the U and V displacement, respectively. B ecause the displacements constantly vary throughout the RVE, the average overall displacements were used for each direction. U sing this method, the displacements that were used for the model in the U and V directi ons were 12.33 m and 21.75 m respectively. T hose values were achieved by taking the average displacement occurring in each direction along the length of the RVE. PAGE 107 107 R esult of R esidual S train and S tress from FE M odel Fig ure 41 0 shows the normal strain for the surface elements of matrix within the RVE calculated in ABAQUS. When comparing the corresponding plots in Fig ure 4 7 the overall trends throughout the surface are in good agreement. The values differ slightly in mo st regions, but significantly for the area along the fiber bundles. It would not be expected for the results to match perfectly as the FE model didnt account for any chemical shrinkage that would be occurring in the actual specimen. Another possible reaso n that could cause such a difference would be the amount of resin in the specimen as opposed to that in the FE model. Having low resin content over the fibers would yield higher strains due to the lack of material to support the surrounding compressive str ains. Therefore, it does seem reasonable that such a large discrepancy would occur in the values for the maximum tensile strains. The displacement information obtained from experiment is inputted into the FE model as the PBCs. The residual stress contours are shown in Fig ure 411. PAGE 108 108 A B Figure 4 1 0 Strain map from FE modeling. A) Strain xx B) Strain yy PAGE 109 109 A B Figure 4 1 1 Residual stress maps from FE modeling. A) Stressxx B) Stressyy PAGE 110 110 From the contours above, it s seen that the residual stress maps have similar color contour maps as residual strains. However, the residual stresses on the surface are all tension with distribution variat ion. P eak tensile values are occurring over the fiber bundles and the smaller values are within the resin rich areas The peak tensile stress was approximately 2 9 MPa and that was located in the very center of the RVE over the fiber bundle where the resin c ontent is lowest. A lthough the values are relatively low, those are the stresses within the resin on the surface. According to the data sheet provided by the prepreg manufacturer, the tensile strength of the resin is 82.7MPa. Although the actual stresses d eveloped were below that value, they are still important to be underst oo d since they limit the usable strength remaining in the structure. Summary The improved CRM was utilized to measure the process induced residual strains for the plain weave composite specimen. An automated strain analysis system was developed to calculate the full field residual strains throughout the entire region of interest. Significant strains were found within the RVE, both in compression and in tension. After examining the strain maps for several specimens, it is seen that fairly repeatable strain values were obtained in the extreme strain regions. A FE model was developed to simulate the curing cycle and the results from the FE model shows similar overall trends as the experiment data. The displacement information obtained from the experiment was inputted to the FE model as PBCs. Residual stresses corresponding to the residual strains were predicted. Although the residual stresses developed were small, they are still important to be underst oo d since they limit the usable strength remaining in the structure. PAGE 111 111 Shrinkage M easurement of C oncrete M aterial CRM was used here to determine the shrinkage that develops in concrete during the curing process. A moir interferometry diffraction grating was replicated on the specimen during the curing process and the specimen was stored in the chamber which controlled specific drying condition for six days after it was demolded. Shrinkage as a function of time, humidity and temperature was measured. During this period, the specimen was removed from the chamber and a quick measurement of shrinkage was made using moir Interferometry every 24 hours. C onsecutive phase shifted fringe patterns were recorded by a digital camera and analyzed using the aut omated strain analysis system to obtain the full field strain information. Introduction Shrinkage is one of the fundamental characteristics of concrete which mainly results from moisture diffusion and self desiccation. Moreover, shrinkage depends on many factors, such as water/cement ratio, surrounding relative humidity and temperature, aggregate quantity, and the size and shape of specimen. This crucial characteristic influences the development of stress and induces subsequent cracking of concrete. In ord er to better assess the structural life of concrete, the shrinkage, which is also time dependent, must be measured in different drying conditions. Much research has been devoted to determining the shrinkage that develops in concrete. ASTM C157 is the standard test method for measuring the length change of concrete using comparator [129] However, this method only measured the averaged shrinkage of concrete and did not provide information on the local area. In 1998, J.K. Kim and C.S. Lee used embedded strain gauges to measure internal shrinkage at the specific location [130] concrete PAGE 112 112 samples using digital photogrammetric method [131] Despite a full field measurement, this optical technique does not have high resolution. The Curing Reference Method (CRM) was invented in 1996 to determine full field residual strains that develop as the laminated composite cools via m oir Interferometry [124] It also could be applied successfully on woven composites [132 134] Moreover, post gel chemical shrinkage of epoxy could be determined in the similar principle [133 135] By using CRM, the in plane deformation of the concrete can be recorded on a full field scale. The use of CRM does not reinforce the specimen and avoid destruction of the specimen. Additionally, the grating does not degrade as time. Therefore, time dependent deformation could be measured. Specimen P reparation To determine the shrinkage, the full field of displacement is obtained via m oir Interferometry and the Curing Reference Method. Chen et al. [136] developed an effective procedure to prepare the specimen and take the measurement in the ESA lab. First, a diffraction grating was replicated onto the specimen from a submaster grating while the specimen was at a free stress state. After the specimen was separated from the submaster grating, the submaster grating was used as the reference grating for adjustment of moir interferometer to obtain null field. Then the specimen was mounted before m oir interferometer and the initial fringe pattern was recorded to determine the initial deformation. Then the specimen was stored in the chamber which provided a specific drying condition for six days. During this period, every 24 hours the specimen was removed from the chamber and a quick measurement was taken to observe how the shrinkage responds to the dryi ng condition. The diffraction grating used had a PAGE 113 113 frequency of 1200 lines/mm. The detail of the procedure for the grating replication can be found in reference [136] In the experiment, a 100% air tight food storage container was used as the chamber. Suffic ient desiccates were put into the container to create nearly 0% relative humidity (RH) inside. On the other hand, 100% RH was easily created by replacing desiccates with some water. 50% RH is room relative humidity. As to temperature control, the oven and Result for C ement w ithout A ggregate The cement paste specimens with w/c ratio=0.5 were used for the testing to determine the shrinkage under different drying environment. For each test, both U and V field were recorded to show similar deformation in both x and y direction because of symmetric condition. Fig ure 4 1 2 ~1 3 show s the fringe PAGE 114 114 Figure 41 2 U field under 50% RH and 23oC Figure 41 3 V field under 50% RH and 23oC Automated strain analysis system was used to analyze the fringe patterns to obtain the full field displacements and strains for each day. Fig ure 4 1 4 ~1 5 only show s the full field displacements on day2 for PAGE 115 115 Displacement: U field (m) Pixel Pixel 0 200 400 600 800 1000 0 200 400 600 800 1000 6 4 2 0 2 4 6 Displacement: V field (m) Pixel Pixel 0 200 400 600 800 1000 1000 800 600 4000 200 0 6 4 2 0 2 4 6 Figure 41 4 Displacement field in the Rectangular Coordinate system for day 2 Displacement: R field (m) PixelPixel 0 200 400 600 800 1000 0 200 400 600 800 1000 6 4 2 0 Displacement: field (m) PixelPixel 0 200 400 600 800 1000 0 200 400 600 800 1000 0.3 0.2 0.1 0 0.1 0.2 Figure 41 5 Displacement field in the Polar Coordinate system for day 2 Fig ure 4 1 6 ~1 7 show s the corresponding strain field in both coordinate systems Again, the results from U and V field are very similar because of the symmetry. PAGE 116 116 Strain xx ( ) Pixel Pixel 0 200 400 600 800 1000 0 200 400 600 800 1000 1200 1100 1000 900 800 700 600 500 400 300 200 Strain yy ( ) Pixel Pixel 0 200 400 600 800 1000 0 200 400 600 800 1000 1200 1100 1000 900 800 700 600 500 400 300 200 Shear Strain xy ( ) Pixel Pixel 0 200 400 600 800 1000 1000 800 600 400 200 0 300 200 100 0 100 200 300 Figure 41 6 Strain field in the Rectangular Coordinate system for day 2 PAGE 117 117 Strain rr ( ) Pixel Pixel 0 200 400 600 800 1000 0 200 400 600 800 1000 1200 1100 1000 900 800 700 600 500 400 300 200 Strain ( ) Pixel Pixel 0 200 400 600 800 1000 0 200 400 600 800 1000 550 500 450 400 350 300 250 200 150 Strain r ( ) Pixel Pixel 0 200 400 600 800 1000 0 200 400 600 800 1000 80 60 40 20 0 20 40 60 Figure 41 7 Strain field in the Polar Coordinate system for day 2 Fig 4 1 8 shows the shrinkage data from the center of the specimen to the expose surface from day 2 to day area close to the exposed surface than that in the center. This is due to faster drying in the area close to the exposed surface. For RH=0%, the result has the same trend as RH=50% except for lar ger shrinkage. PAGE 118 118 Figure 41 8 Shrinkage from center to expose surface In the same experimental and analytical technique, temperature effect was also investigated. The fringe patterns are not shown here. Fig 4 19~ 2 0 shows the results for temperature and relative humidity effect on the shrinkage of concrete specimens in average sense. Figure 419. Temperature effect PAGE 119 119 Figure 42 0 Relative humidity effect Result for C ement with C oarse A ggregate Gravel was embedded into cement paste (w/c=0.5) clos e to the grating, but did not contact with the grating. The arrangement of the gravel is shown as Fig ure 4 2 1 with a side view and top view. The specimen was stored in relative humidity 0% and room temperature. Silicone rubber mold Concrete Epoxy Silicone rubber grating Gravel Figure 4 2 1 The side and top views of gr avel test Moir fringe patterns from Day 1 to Day 5 were recorded. O nly day 5 moir fringe patterns are shown here in Fig ure 4 2 2 There were some areas of lower fringe density marked by blue boxes. These areas corresponded to the locations of the embedded PAGE 120 120 gravels. A crack occurred within the area of the red circle. This was due to the constraint from the embedded gravel. Fig ure 42 2 Day 5 moir fringe patterns for the gravel test Normal strain analysis was performed in Fig ure 4 2 3 Clearly, the locati ons of gravel show lower shrinkage. The plots just displayed the normal strain distributions along the horizontal centerline in U field and the vertical centerline in V field. More result for cement with fine aggregate can be found in [137] Strain xx ( ) Pixel Pixel 0 200 400 600 800 1000 1000 800 600 400 200 0 3000 2500 2000 1500 1000 0 100 200 300 400 500 600 700 800 900 1000 3500 3000 2500 2000 1500 1000 500 Position Shrinkage Strain yy ( ) 0 200 400 600 800 1000 1000 800 600 400 200 0 3000 2500 2000 1500 1000 0 200 400 600 800 1000 3500 3000 2500 2000 1500 1000 500 Position Shrinkage Fig ure 42 3 The normal s train analysis for gravel test PAGE 121 121 Summary The experimental technique for measuring the shrinkage of concrete has been developed based on curing reference method and moir interferometry. The automated strain analysis system was applied to obtain the full field displacement and strain maps successfully T he experimental data show that t he shrinkage of concrete increases with time and is very sensitive to surrounding environment. I t increases with t emperature. A nd a decrease in rel ative humidity also increases the shrinkage. This experimental method was also employed to investigate the effect of coarse or fine aggregate on shrinkage of concrete. Different quantities of sands can be added to cement paste to see how shrinkage behaves. A numerical model can be developed in order to obtain material properties from the complex geometry used in our test. Material P roperty I dentification of L aminate U sing F ull F ield M easurements This section covers the cooperation work with Chiristian Gogu. The objective is to identify the four ply elastic constants 121212,,, EEG from the moir full field displacement measurements on a laminate plate with a hole using Bayesian identification approach. Only the experiment part is covered here. All the details of material property identification can be found in [138] Intro duction Identification of the four orthotropic elastic constants of a composite from a tensile test on a plate with a hole was carried out by several authors in the past [139 140] using finite element model updating based on a least squares framework. Some of the advantages of doing the identification from measurements on a plate with a hole are the ability to identify all four properties at the same time from a single experiment. This is possible because of the heterogeneous strain field exhibited during this test which PAGE 122 122 involves all four elastic constants, unlike tension tests on simple rectangular specimen which exhibit uniform strain fields and usually involve only two of the elastic constants. A further advantage is that the heterogeneous fields can provide information on spatial variations of the material properties as well. Moir interferometry was used to capture the displacement nonuniformity on the specimen used for identification because of its high spatial resolution. The objective is to use Bayesian identification to for identif ic ation of the four ply elastic constants from the fullfield displacement measurements on a laminate plate with a hole. The interest of identifying ply properties is that it allows to obtain both the extensional and the bending stiffnesses of the laminate. Due however to the varying sensitivity of the strain and displacement fields to the different ply properties, it is of primary importance to estimate the uncertainty with which these properties are identified. The Bayesian approach developed by Gogu [138] will allow us to do this by taking into account the physics of the problem (i.e. the different sensitivities of the fields to the different properties), measurement uncertainty as well as uncertainty on other input parameters to the model such as the specimen geometry. The experiment considered for the identification is a tensile test on a composite plate with a hole. The laminate is made of graphite/epoxy with a stacking sequence of 45,45,0s The plate has the dimensions given in Figure 42 4 with a ply thickness of 0.16 mm. The applied tensile force is 1200 N. The fullfield measurement takes place on a square area 20 x 20 mm2 around center of the hole. PAGE 123 123 Figure 42 4 Specimen geometry. The specimen material is graphite/epoxy and the stacking sequence 45,45,0s Experiment T he dimension of the tes t specimen is slightly different from the designed. The width of the specimen was 24.3 mm, the diameter of the hole was 4 10 mm and the total laminate thickness was 0.78mm. The plate was made out of a Toray T800/3631 graphite/epoxy prepreg. The manufacturers specifications for this material are given in Table 4 2 together with the properties obtained by Noh [141] Noh obtained the material properties based on a four points bending test at the University of Florida on a composite made from the exact same prepreg roll that we used. Table 42. Manufacturers specifications and properties found by Noh (2004) based on a four points bending test Para meter E 1 (GPa) E 2 (GPa) G 12 (GPa) v 12 Manufactures specification 162 7.58 4.41 0.34 Noh (2004) values 144 7.99 7.78 0.34 A picture of our specimen with the diffraction grating ( 1200 lines/mm) is shown in Figure 4 2 5 The specimen was loaded at 700 N for the measurements to consider the safety of the material. PAGE 124 124 Figure 42 5 Test specimen with grating (1200 lines/mm) The experiment setup is shown in Figure 42 6 PEMI II 2020 X Moir interferometer using a Pulnix TM 1040 digital camera were utilized to ca pture phase shifted fringe patterns. MTI 30K machine was used to apply tension load. Rotations of the grips holding the specimen were allowed by using a lubricated ball bearing for the bottom grip and two lubricated shafts for the top grip. This allowed reducing parasitic bending during the tension test. Figure 42 6 Experimental setup PAGE 125 125 Fringe P atterns and R esult The fringe patterns observed for a force of 700 N are shown in Figures 4 2 7 The two smaller holes in the fringe patterns and other parasitic lines are due to imperfections in the diffraction grating such as air bubbles and dust from the cover magic Scott tapes. A B Figure 4 2 7 Fringe patterns A) U field B) V field The automated strain analysis system was used to analyze the fringe patter ns and extract the displacement fields. The corresponding displacement maps are shown in Figures 42 8 Both of the U and V field displacement maps will be used in the Bayesian identification procedure. PAGE 126 126 Displacement: U field (m) 100 200 300 400 500 600 700 100 200 300 400 500 600 700 6 4 2 0 2 4 A Displacement: V field (m) 100 200 300 400 500 600 700 100 200 300 400 500 600 700 6 4 2 0 2 4 6 B Figure 42 8 Displacement maps A) U f ield B) V field The results of the identification are provided in Tables 4 3 Table 43. Mean values and coefficient of variation of the identified posterior distribution from the Moir interferometry experiment Parameter E 1 (GPa) E 2 (GPa) G 12 (GPa) v 12 Mea n value 138 7.48 5.02 0.33 COV (%) 3.12 9.46 4.29 10.3 Although strain maps are not necessary for the material property identification procedure, the strain maps obtained from the automated strain analysis system are shown in Figure 429. They have good agreement with those obtained in FE simulation for most the areas except the extremely large values close to the edges. PAGE 127 127 Strain xx ( ) 100 200 300 400 500 600 700 100 200 300 400 500 600 700 1600 1400 1200 1000 800 600 400 200 0 200 400 A Strain yy ( ) 100 200 300 400 500 600 700 100 200 300 400 500 600 700 0 200 400 600 800 1000 B Figure 429. Strain maps A) Strain xx B) Strain yy Summary Overall th e mean values of the identified distribution are in agreement with the manufacturers specifications. The largest difference is in longitudinal Youngs modulus, which could seem somewhat surprising. However Noh[141] found a similar value on the exact same prepreg roll that we used. The mean values of E212 and G12 are close to the specification values. G12 is far however from Nohs values but it should be noted that the four point bending test is relatively poor for identifying G12. S ummary I n the present chapter, selected three applications of Auto mated Strain Analysis System in the ESA lab were introduced and the results were illustrated. T he system together with moir interferometry can give accurate measurement of residual strain of textile composite and shrinkage of concrete material and provide accurate displacement data for material property identification of laminate. PAGE 128 128 CHAPTER 5 CONCLUSIONS AND SUGGESTED FUTURE WORK Conclusions Moir interferometry has been widely applied for industrial and scientific studies because of its capability of f ullfield displacement measurement, high sensitivity and high spatial resolution. T he output of m oir is represented as interference fringe patterns which need to be analyzed to obtain the desired physical parameters such as displacement, strain, etc. Trad itionally, the fringe pattern was analyzed manually via drawing the gage lines on some points and counting the number of fringes crossing them To improve the efficiency of the experiments using moir interferometry and enhance the accuracy of the fringe p attern analysis, an automated strain analysis system was developed in this dissertation. The major contributions in this dissertation include: Complete survey and investigation on the existing image processing techniques for fringe pattern analysis; appro priate algorithms were adopted and improved. Implementation of the quality guided phase unwrapping Improvement of the quality guided phase unwrapping; enhancement of the efficiency of the algorithm for application of large image Development of phase repair for the inconsistent areas Development of automatic detection of the inconsistent areas Development of the hybrid method O/DFM. Development of global surface fit to calculate strain effectively Development of local surface fit to calculate strain effectively PAGE 129 129 Development of a Windows GUI based software Automated Strain Analysis System using MatLab 2009 for automatic fringe pattern analysis. This software includes all the algorithms covered in the dissertation. It also includes some useful tools for post processing of the data. Application of the system for residual stress characterization of plain woven textile composite Application of the system for shrinkage measurement of concrete materials without and with aggregate. Application of the system f or material property identification of laminate using openhole test. Suggest ed F uture W ork S o far the automated strain analysis system has been successfully applied for several projects in the ESA lab. A lthough some preliminary results were obtained, the following tasks may be completed in the future, Enhancement of the computational efficiency of global surface fit using thin plate spline. With the increased number of seed points on those areas close to the edge, it will calculate the strain within thos e areas better. Improvement of the algorithm for the automatic detection of inconsistent areas in the unwrapped phase map PAGE 130 130 APPENDIX A BILINEAR FUNCTION FO R LOCAL SURFACE FIT BASED STRAIN CALCULA TION Bilinear function can be used to fit the 2D data withi n small patches. The unwrapped phase is assumed to have the bilinear distribution within one patch as shown in Equation (A 1) 1234, fijccicjcij ( A 1) where fij is the smoothed unwrapped phase. ij is the image coordinate system with unit of pixel. 1234,,, cccc are the coefficients calculated from the surface fit. All the points within the patch of size mn pixels are used to for m the equation group, 11112131411 22212232422 1234, ...... ,mnmnmn mnmnmnmnfijccicjcij fijccicjcij fijccicjcij ( A 2) ,4 mnB 4,1C and ,1 mnY matrix can be formed, 1111 2222 ,41 1 ...... 1mn mnmnmnmnijij ijij B ijij ( A 3) 1 2 4,1 3 4c c C c c ( A 4) PAGE 131 131 111 222 ,1, ...... ,mn mnmnmnfij fij Y fij ( A 5) Then the coefficients can be calculated as, 1 TT CBBBY ( A 6) Finally the gradients (strains) can calculated analytically as, 34 24, ff ijcci xj ff ijccj yi ( A 7) The negative sign expression of / fy comes from the opposite direction in image system and rectangular coordinate system. PAGE 132 132 APPENDIX B CUBIC FUNCTION FOR L OCAL SURFACE FIT BAS ED STRAIN CALCULATIO N Cubic function can be used to fit the 2D data within small patches. The unwrapped phase is assumed to have the cubic distribution within one patch as shown in Equation (B 1) 223223 12345678910, fijccicjcijcicjcicijcijcj ( B 1 ) where fij is the smoothed unwrapped phase. ij is the image coordinate system with unit of pixel. 12345678910,,,,,,,,, cccccccccc are the coefficients calculated from the surface fit. All the points within the patch of size mn pixels are used to form the equation group, 223223 11112131411516171811911101 223223 22212232422526272822922102 223 123456 8, ...... ,mnmnmn mnmnmnmnmnmnmnfijccicjcijcicjcicijcijcj fijccicjcijcicjcicijcijcj fijccicjcijcicjcici 223 9 10 mnmnmnmnmnjcijcj ( B 2 ) ,10 mnB 10,1C and ,1 mnY matrix can be formed, 223223 111111111111 223223 222222222222 ,10 2232 231 1 ...... 1mn mnmnmnmnmnmnmnmnmnmnmnmnijijijiijijj ijijijiijijj B ijijijiijijj ( B 3 ) PAGE 133 133 1 2 3 4 5 10,1 6 7 8 9 10c c c c c C c c c c c ( B 4 ) 111 222 ,1, ...... ,mn mnmnmnfij fij Y fij ( B 5 ) Then the coefficients can be calculated as, 1 TT CBBBY ( B 6 ) Finally the gradients (strains) can calculated analytically as, 22 3468910 22 245789,223 232 ff ijccicjcicijcj xj ff ijccjcicicijcj yi ( B 7 ) The negative sign expression of / fy comes from the opposite direction in image system and rectangular coordinate system. PAGE 134 134 LIST OF REFERENCES 1. D. Post, B. Han and P. Ifju, High Sensitivity Moire: Experimental Analysis for Mechanics and Materials Springer Verlag, New York (1994). 2. A. Livnat and D. Post, "The Governing Equations for Moire Interferometry and Their Identity to Equations of Geometrical Moire," Exp Mech 25(4), 360366 (1985) 3 G. Peacock, Miscellaneous Works of the Late Thomas Young John Murray Publishers, London (1855). 4. J. Dyson, "The Rapid Measurement of Photographic Records of Interference Fringes," Appl Optics 2(5), 487489 (1963) 5. G. D. Dew, "Method for Precise Eval uation of Interferograms," J Sci Instrum 41(3), 1601 65 (1964) 6. A. F. Bastawros and A. S. Voloshin, "Thermal Strain Measurements in Electronic Packages through Fractional Fringe Moire Interferometry," Journal of Electronic Packaging, Trans. ASME 112(4), 30 3 308 (1990) 7. T. Kreis, holographic interferometry Academic Verlag, Berlin (1996). 8. D. Malacara, M. Servin and Z. Malacara, Interferogram analysis for optical testing Marcel Dekker Inc., New York (1998). 9. A. Jain, Fundamentals of digital image proc essing, PrenticeHall International Inc., London (1989). 10. G. Jin, Computer aided optical measurements Tsinghua University Press, Beijing (1997). 11. T. Yatagai, "Intensity based analysis methods. In Interferogram Analysis (Eds D.W. Robinson and G.T. Reid)," 7393 (1993) 12. T. Yatagai, S. Nakadate, M. Idesawa and H. Saito, "Automatic Fringe Analysis Using Digital Image Processing Techniques," Opt Eng 21(3), 432435 (1982) 13. G. A. Mastin and D. C. Ghiglia, "Digital Extraction of Interference Fringe Con tours," Appl Optics 24(12), 17271728 (1985) 14. T. Y. Chen and C. E. Taylor, "Computerized Fringe Analysis in Photomechanics," Exp Mech 29(3), 323329 (1989) 15. F. Bertolino and F. Ginesu, "Semiautomatic Approach of Grating Techniques," Opt Laser Eng 19( 4 5), 313323 (1993) 16. A. E. Ennos, D. W. Robinson and D. C. Williams, "Automatic Fringe Analysis in Holographic Interferometry," Opt Acta 32(2), 135145 (1985) PAGE 135 135 17. S. Krishnaswamy, "Algorithm for Computer Tracing of Interference Fringes," Appl Optics 30 (13), 16241628 (1991) 18. J. A. Aparicio, J. L. Molpeceres, A. M. Defrutos, C. Decastro, S. Caceres and F. A. Frechoso, "Improved Algorithm for the Analysis of Holographic Interferograms," Opt Eng 32(5), 963969 (1993) 19. D. S. Zhang, M. Ma and D. D. Aro la, "Fringe skeletonizing using an improved derivative sign binary method," Opt Laser Eng 37(1), 5162 (2002) 20. D. S. Zhang, D. D. Arola and M. Ma, "Digital image processing of moire fringe patterns with an application to fractures in bovine dentin," Opt Eng 41(5), 11151121 (2002) 21. B. T. Han, "Interferometric Methods with Enhanced Sensitivity by Optical Digital Fringe Multiplication," Appl Optics 32(25), 47134718 (1993) 22. A. Rosenfeld, "Characterization of Parallel Thinning Algorithms," Inform Cont rol 29(3), 286291 (1975) 23. A. Rosenfeld and A. Kak, Digital picture processing, Vol. I Academic Press, New York (1982). 24. J. B. Schemm and C. M. Vest, "Fringe Pattern Recognition and Interpolation Using Non Linear RegressionAnalysis," Appl Optics 22 (18), 28502853 (1983) 25. A. Colin and W. Osten, "Automatic Support for Consistent Labeling of Skeletonized Fringe Patterns," J Mod Optic 42(5), 945954 (1995) 26. W. W. Macy, "Two Dimensional FringePattern Analysis," Appl Optics 22(23), 38983901 (1983) 27. P. L. Ransom and J. V. Kokal, "Interferogram Analysis by a Modified Sinusoid Fitting Technique," Appl Optics 25(22), 41994204 (1986) 28. T. Yatagai, "Automated Fringe Analysis Techniques in Japan," Opt Laser Eng 15(2), 7991 (1991) 29. M. Takeda, H. Ina and S. Kobayashi, "Fourier Transform Method of FringePattern Analysis for Computer Based Topography and Interferometry," J Opt Soc Am 72(1), 156160 (1982) 30. K. H. Womack, "Frequency Domain Description of Interferogram Analysis," Opt Eng 23(4), 396400 (1984) 31. D. J. Bone, H. A. Bachor and R. J. Sandeman, "FringePattern Analysis Using a 2 D Fourier Transform," Appl Optics 25(10), 16531660 (1986) PAGE 136 136 32. T. Kreis, "Digital Holographic Interference Phase Measurement Using the Fourier Transform Method," J Opt Soc Am A 3(6), 847855 (1986) 33. C. Roddier and F. Roddier, "Interferogram Analysis Using Fourier Transform Techniques," Appl Optics 26(9), 16681673 (1987) 34. R. J. Green, J. G. Walker and D. W. Robinson, "Investigation of the Fourier Transform Method of Fringe PatternAnalysis," Opt Laser Eng 8(1), 2944 (1988) 35. K. Freischlad and C. L. Koliopoulos, "Fourier Description of Digital Phase Measuring Interferometry," J Opt Soc Am A 7(4), 542551 (1990) 36. D. Tentori and D. Salazar, "Hologram Inter ferometry Carrier Fringes," Appl Optics 30(35), 51575158 (1991) 37. P. J. Bryanstoncross, C. Quan and T. R. Judge, "Application of the Fft Method for the Quantitative Extraction of Information from HighResolution Interferometric and Photoelastic Data," Opt Laser Technol 26(3), 147155 (1994) 38. J. Gu and F. Chen, "Fast Fourier Transform, Iteration, and Least Squares Fit Demodulation ImageProcessing for Analysis of SingleCarrier Fringe Pattern," J Opt Soc Am A 12(10), 21592164 (1995) 39. S. De Nicola and P. Ferraro, "Fourier transform method of fringe analysis for moire interferometry," J Opt a Pure Appl Op 2(3), 228233 (2000) 40. J. H. Massig and J. Heppner, "Fringepattern analysis with high accuracy by use of the Fourier transform method: theory and experimental tests," Appl Optics 40(13), 20812088 (2001) 41. L. L. Deck, "Fourier transform phase shifting interferometry," Appl Optics 42(13), 23542365 (2003) 42. O. Y. Kwon, "Multichannel PhaseShifted Interferometer," Opt Lett 9(2), 5961 (1984) 43. B. Brueckmann and W. Thieme, "Computer Aided Analysis of Holographic Interferograms Using the Phase Shift Method," Appl Optics 24(14), 21452149 (1985) 44. M. Chang, C. P. Hu, P. Lam and J. C. Wyant, "HighPrecision Deformation Measurement by Digital Phase Shifting Holographic Interferometry," Appl Optics 24(22), 37803783 (1985) 45. Y. Ishii, J. Chen and K. Murata, "Digital Phase Measuring Interferometry with a Tunable Laser Diode," Opt Lett 12(4), 233235 (1987) 46. M. Kujawinska, "Use of PhaseStepping Automatic Fringe Analysis in Moire Interferometry," Appl Optics 26(22), 47124714 (1987) PAGE 137 137 47. K. Kinnstaetter, A. W. Lohmann, J. Schwider and N. Streibl, "Accuracy of Phase Shifting Interferometry," Appl Optics 27(24), 50825089 (1988) 48. E. Vikhagen and O. J. Lokberg, "Detection of Defects in CompositeMaterials by Television Holography and ImageProcessing," Mater Eval 48(2), 244248 (1990) 49. G. Lai and T. Yatagai, "Generalized Phase Shifting Interferometry," J Opt Soc Am A 8(5), 822827 (1991) 50. K. Okada, A. Sato and J. Tsujiuchi, "Simultaneous Calculation of Phase Distribution and Scanning Phase Shift in Phase Shifting Interferometry," Opt Commun 84(34), 118124 (1991) 51. C. Y. Poon, M. Kujawinska and C. Ruiz, "Automated Fringe Pattern Analysis fo r Moire Interferometry," Exp Mech 33(3), 234241 (1993) 52. C. Y. Poon, M. Kujawinska and C. Ruiz, "Spatial Carrier Phase Shifting Method of Fringe Analysis for Moire Interferometry," J Strain Anal Eng 28(2), 7988 (1993) 53. Y. Morimoto and M. Fujisawa, Fringe Pattern Analysis by a Phase Shifting Method Using Fourier Transform," Opt Eng 33(11), 37093714 (1994) 54. Y. Surrel, "Design of algorithms for phase measurements by the use of phase stepping," Appl Optics 35(1), 5160 (1996) 55. G. Stoilov and T. D ragostinov, "Phase stepping interferometry: Fiveframe algorithm with an arbitrary step," Opt Laser Eng 28(1), 61 69 (1997) 56. G. H. Kaufmann, "Phase shifting in speckle photography," Opt Laser Eng 29(23), 117123 (1998) 57. A. Davila, P. D. Ruiz, G. H. Kaufmann and J. M. Huntley, "Measurement of subsurface delaminations in carbon fibre composites using highspeed phase shifted speckle interferometry and temporal phase unwrapping," Opt Laser Eng 40(56), 447458 (2003) 58. H. M. Xie, Z. W. Liu, D. N. Fan g, F. L. Dai, H. J. Gao and Y. P. Zhao, "A study on the digital nanomoire method and its phase shifting technique," Meas Sci Technol 15(9), 17161721 (2004) 59. H. M. Xie, H. X. Shang, F. L. Dai, B. Li and Y. M. Xing, "Phase shifting SEM moire method," Op t Laser Technol 36(4), 291 297 (2004) 60. M. C. Brito, J. M. Alves, J. M. Serra, R. M. Gamboa, C. Pinto and A. M. Vallera, "Measurement of residual stress in EFG ribbons using a phaseshifting IR photoelastic method," Sol Energ Mat Sol C 87(14), 311316 ( 2005) PAGE 138 138 61. Y. Morimoto, T. Nomura, M. Fujigaki, S. Yoneyama and I. Takahashi, "Deformation measurement by phase shifting digital holography," Exp Mech 45(1), 6570 (2005) 62. J. L. Chen, C. R. Sun, Y. W. Qin and X. H. Ji, "Damage evaluation of subsurface de fect in sandwich by phaseshifting digital shearography," Fracture and Strength of Solids Vi, Pts 1 and 2 306308(399404 (2006) 63. F. Henault, "Conceptual design of a phaseshifting telescopeinterferometer," Opt Commun 261(1), 3442 (2006) 64. Y. Min, M Hong, Z. Xi and L. Jian, "Determination of residual stress by use of phase shifting moire interferometry and holedrilling method," Opt Laser Eng 44(1), 68 79 (2006) 65. C. Quan, S. H. Wang and C. J. Tay, "Nanoscale surface deformation inspection using F FT and phase shifting combined interferometry," Precis Eng 30(1), 2331 (2006) 66. H. N. Yen, D. M. Tsai and J. Y. Yang, "Full field 3 D measurement of solder pastes using LCD based phase shifting techniques," Ieee T Electron Pack 29(1), 5057 (2006) 67. C Han and B. Han, "Phase shifting in achromatic moir interferometry system," OPTICS EXPRESS 15(16), 99709976 (2007) 68. J. M. Huntley, "Automated fringe pattern analysis in experimental mechanics: a review," J Strain Anal Eng 33(2), 105125 (1998) 69. Z. Y. Wang, "Development and Application of Computer aided Fringe Analysis," in Mechanical Engineering, University of Maryland, Maryland (2003). 70. K. Itoh, "Analysis of the Phase Unwrapping Algorithm," Appl Optics 21(14), 24702470 (1982) 71. R. M. Goldste in, H. A. Zebker and C. L. Werner, "Satellite Radar Interferometry Two Dimensional Phase Unwrapping," Radio Science 23(4), 713720 (1988) 72. J. M. Huntley, "NoiseImmune Phase Unwrapping Algorithm," Appl Optics 28(16), 32683270 (1989) 73. D. J. Bone, Fourier Fringe Analysis the 2Dimensional Phase Unwrapping Problem," Appl Optics 30(25), 36273632 (1991) 74. N. H. Ching, D. Rosenfeld and M. Braun, "Twodimensional phase unwrapping using a minimum spanning tree algorithm," Ieee T Image Process 1(3), 3 55 365 (1992) 75. T. R. Judge and P. J. Bryanstoncross, "A Review of Phase Unwrapping Techniques in Fringe Analysis," Opt Laser Eng 21(4), 199239 (1994) PAGE 139 139 76. D. C. Ghiglia and L. A. Romero, "Robust 2Dimensional Weighted and Unweighted Phase Unwrapping That Uses Fast Transforms and Iterative Methods," J Opt Soc Am A 11(1), 107117 (1994) 77. J. A. Quiroga, A. Gonzalezcano and E. Bernabeu, "PhaseUnwrapping Algorithm Based on an Adaptive Criterion," Appl Optics 34(14), 25602563 (1995) 78. D. Kerr, G. H. Kau fmann and G. E. Galizzi, "Unwrapping of interferometric phasefringe maps by the discrete cosine transform," Appl Optics 35(5), 810816 (1996) 79. D. C. Ghiglia and L. A. Romero, "Minimum L(p) norm twodimensional phase unwrapping," J Opt Soc Am A 13(10), 19992013 (1996) 80. T. J. Flynn, "Two dimensional phase unwrapping with minimum weighted discontinuity," J Opt Soc Am A 14(10), 26922701 (1997) 81. D. C. Ghiglia and M. D. Pritt, Two Dimensional Phase Unwrapping: Theory, Algorithms, and Software, Wiley, New York (1998). 82. G. H. Kaufmann, G. E. Galizzi and P. D. Ruiz, "Evaluation of a preconditioned conjugategradient algorithm for weighted least squares unwrapping of digital specklepattern interferometry phase maps," Appl Optics 37(14), 30763084 (1998) 83. W. Xu and I. Cumming, "A regiongrowing algorithm for InSAR phase unwrapping," Ieee Transactions on Geoscience and Remote Sensing 37(1), 124134 (1999) 84. M. A. Herraez, D. R. Burton, M. J. Lalor and M. A. Gdeisat, "Fast two dimensional phase unwrap ping algorithm based on sorting by reliability following a noncontinuous path," Appl Optics 41(35), 74377444 (2002) 85. W. S. Li and X. Y. Su, "Phase unwrapping algorithm based on phase fitting reliability in structured light projection," Opt Eng 41(6), 1 3651372 (2002) 86. J. Arines, "Least squares modal estimation of wrapped phases: application to phase unwrapping," Appl Optics 42(17), 33733378 (2003) 87. A. Baldi, "Phase unwrapping by region growing," Appl Optics 42(14), 24982505 (2003) 88. X. Y. Su a nd W. J. Chen, "Reliability guided phase unwrapping algorithm: a review," Opt Laser Eng 42(3), 245261 (2004) 89. J. M. Bioucas Dias and G. Valadao, "Phase unwrapping via graph cuts," Lect Notes Comput Sc 3522 360 367 (2005) 90. Y. G. Lu, X. Z. Wang and G. T. He, "Phase unwrapping based on branch cut placing and reliability ordering," Opt Eng 44(5), (2005) PAGE 140 140 91. B. Marendic, Y. Y. Yang and H. Stark, "Phase unwrapping using an extrapolationprojection algorithm," J Opt Soc Am A 23(8), 18461855 (2006) 92. Y. J. Zhu, L. R. Liu, Z. Luan and A. H. Li, "A reliable phase unwrapping algorithm based on the local fitting plane and quality map," J Opt a Pure Appl Op 8(6), 518523 (2006) 93. Y. G. Lu, X. Z. Wang and X. P. Zhang, "Weighted least squares phase unwrapping algorithm based on derivative variance correlation map," Optik 118(2), 6266 (2007) 94. L. An, Q. S. Xiang and S. Chavez, "A fast implementation of the minimum spanning tree method for phase unwrapping," Ieee T Med Imaging 19(8), 805808 (2000) 95. W. Bos saert, R. Dechaene and A. Vinckier, "Computation of finite strains from moire displacement patterns," J Strain Anal Eng 3 63 75 (1968) 96. H. A. Vrooman and A. A. M. Maas, "ImageProcessing Algorithms for the Analysis of Phase Shifted Speckle Interference Patterns," Appl Optics 30(13), 16361641 (1991) 97. A. J. Moore and J. R. Tyrer, "Twodimensional strain measurement with ESPI," Opt Laser Eng 24(56), 381402 (1996) 98. J. Morton, D. Post, B. Han and M. Y. Tsai, "A Localized Hybrid Method of Stress Anal ysis a Combination of Moire Interferometry and Fem," Exp Mech 30(2), 195200 (1990) 99. M. Y. Tsai and J. Morton, "New Developments in the Localized Hybrid Method of Stress Analysis," Exp Mech 31 298 305 (1990) 100. J. Rhee, S. He and R. E. Rowlands, "Hy brid Moire numerical stress analysis around cutouts in loaded composites," Exp Mech 36(4), 379387 (1996) 101. K. R. Castleman, Digital Image processing, PrenticeHall, New Jersey (1996). 102. R. C. Gonzalez and R. E. Woods, Digital Image processing Pears on Prentice Hall, New York (2007). 103. M. L. Speriatu, "Temperature Dependent Mechanical Properties of Composite Materials and Uncertainties in Experimental Measurement," in Mechanical and Aerospace Engineering, University of Florida, Gainesville (2005). 104. K. G. Larkin and B. F. Oreb, "Design and Assessment of Symmetrical PhaseShifting Algorithms," J Opt Soc Am A 9(10), 17401748 (1992) PAGE 141 141 105. J. R. Lee, J. Molimard, A. Vautrin and Y. Surrel, "Digital phase shifting grating shearography for experimental analysis of fabric composites under tension," Compos Part a Appl S 35(78), 849859 (2004) 106. J. Duchon, "Splines minimizing rotationinvariant semi norms in Sobolev spaces. In Constructive Theory of Functions of Several Variables. vol 1, pp 85100 (Spri nger Berlin / Heidelberg)," Constructive Theory of Functions of Several Variables 1 85100 (1976) 107. G. Donato and S. Belongie, "Approximate thin plate spline mappings," Computer Vision Eccv 2002 Pt Iii 2352 2131 (2002) 108. I. Barrodale, D. Skea, M. Berkley, R. Kuwahara and R. Poeckert, "Warping Digital Images Using Thin Plate Splines," Pattern Recogn 26(2), 375 376 (1993) 109. H. Guo, J. Y. Jiang and L. M. Zhang, "Building a 3D morphable face model by using thin plate splines for face reconstruction, Advances in Biometric Person Authentication, Proceedings 3338, 258 267 (2004) 110. M. J. D. Powell, "A thin plate spline method for mapping curves into curves in two dimensions," Computational Techniques and Applications (CTAC95), Melbourne, Australia, 1995 (1995) 111. F. Girosi, M. Jones and T. Poggio, "Regularization Theory and Neural Networks Architectures," Neural Comput 7(2), 219269 (1995) 112. R. C. Gonzalez, R. E. Woods and S. L. Eddins, Digital Image processing, Pearson Prentice Hall, New York (2002). 113. R. C. Gonzalez, R. E. Woods and S. L. Eddins, Digital Image processing using MATLAB Pearson Prentice Hall, New York (2004). 114. H. T. Hahn and N. J. Pagano, "Curing Stresses in Composite Laminates," J Compos Mater 9(Jan), 91106 (1975) 115. Z. W u, J. A. Lu and B. T. Han, "Study of residual stress distribution by a combined method of moire interferometry and incremental hole drilling, part I: Theory," J Appl MechT Asme 65(4), 837843 (1998) 116. O. Sicot, X. L. Gong, A. Cherouat and J. Lu, "Deter mination of residual stress in composite laminates using the incremental holedrilling method," J Compos Mater 37(9), 831844 (2003) 117. H. T. Hahn, "Residual Stresses in Polymer Matrix Composite Laminates," J Compos Mater 10(Oct), 266 278 (1976) 118. R. Y. Kim and H. T. Hahn, "Effect of Curing Stresses on the First Ply Failure in Composite Laminates," J Compos Mater 13(Jan), 216 (1979) PAGE 142 142 119. B. Benedikt, R. Rupnowski, L. Kumosa, J. K. Sutter, P. K. Predecki and M. Kumosa, "Determination of interlaminar residual thermal stresses in a woven 8HS graphite/PMR 15 composite using X ray diffraction measurements," Mech Adv Matl Struct 9(4), 375394 (2002) 120. I. G. Zewi, I. M. Daniel and J. T. Gotro, "Residual Stresses and Warpage in Woven Glass Epoxy Laminates," Exp Mech 27(1), 4450 (1987) 121. W. A. Schulz, "Determination of residual stress and thermal behavior for composite laminates," in Mechanical and Aerospace Engineering, University of Florida, Gainesville (2005). 122. X. G. Huang, J. W. Gillespie and T. Bogetti, "Process induced stress for woven fabric thick section composite structures," Compos Struct 49(3), 303 312 (2000) 123. R. Karkkainen, "Stress Gradient Failure Theory for Textile Structural Composites," in Mechanical and Aerospace Engineering, Univ ersity of Florida, Gainesville (2006). 124. P. G. Ifju, X. Niu, B. C. Kilday, S. C. Liu and S. M. Ettinger, "Residual strain measurement in composites using the curereferencing method," Exp Mech 40(1), 2230 (2000) 125. P. G. Ifju, J. E. Masters and W. C. Jackson, "The Use of Moire Interferometry as an Aid to Standard Test Method Development for Textile CompositeMaterials," Compos Sci Technol 53(2), 155163 (1995) 126. N. M. Strickland, "Residual Strain Measurement of Plain Weave Composites Using The Cure Reference Method," in Mechanical and Aerospace Engineering, University of Florida, Gainesville (2007). 127. V. Carvelli and C. Poggi, "A homogenization procedure for the numerical analysis of woven fabric composites," Compos Part aAppl S 32(10), 1425143 2 (2001) 128. J. B. Berman and S. R. White, "Theoretical modelling of residual and transformational stresses in SMA composites," Smart Mater Struct 5(6), 731743 (1996) 129. ASTM, "Standard Test Method for Length Change of Hardened Hydraulic Cement Mortar and Concrete," in C157/C157M 06. 130. J. K. Kim and C. S. Lee, "Prediction of differential drying shrinkage in concrete," Cement Concrete Res 28(7), 985994 (1998) 131. F. Yilmazturk, Kulur, S., Pekmezci, B.Y., "Measurement of Shrinkage in Concrete Samples using Digital Photogrammetric Methods," The International Archives of the Photogrammtry, Remote Sensing and Spatial Information Sciences 34( 2003) PAGE 143 143 132. N. M. Strickland, Yin, W.Q., Ifju, P.G., "Residual Strain Measurement Analysis of Woven Composites usin g Phase Shifting," in SEM Annual Conference and Exposition 10th International Congress and Exhibition on Experimental and Applied Mechanics pp. 18551860 (2007). 133. W. Q. Yin, Strickland, N. M., Ifju, P.G., "Residual Strain Measurement and Analysis of Plain Woven Composites," in 49th AIAA/ASME/ASCE/AHS/ASC Structures, Structural Dynamics, and Materials Conference Chicago (2008). 134. W. Q. Yin, Strickland, N. M., Ifju, P.G., "Residual Strain Measurement and Residual Stress Prediction of Woven Composit es," in Society for Experimental Mechanics 11th International Congress and Exhibition on Experimental and Applied Mechanics pp. 687694, Orlando, FL (2008). 135. S. C. Liu, "Residual Stress Characterization for Laminated Composites," in Mechanical Engineering, University of Florida, Gainesville (1999). 136. T. C. Chen, Yin, W.Q., Ifju, P.G., "Shrinkage Measurement of Concrete Using Phase Shifting," in Society for Experimental Mechanics 11th International Congress and Exhibition on Experimental and Appl ied Mechanics pp. 10921100 (2008). 137. T. C. Chen, Yin, W.Q., Ifju, P.G., "Experimental Investigation of Aggregate Effects on Shrinkage Behavior in Concrete Materials," in Society for Experimental Mechanics 12th International Congress and Exhibition o n Experimental and Applied Mechanics (2009). 138. C. Gogu, "Facilitating Bayesian Identification of Elastic Constants Through Dimensionality Reduction and Response Surface Methodology," in Mechanical and Aerospace Engineering, p. 229, University of Florida, Gainesville (2009). 139. D. Lecompte, Sol, H., Vantomme, J., and Habraken, A. M., "Identification of elastic orthotropic material parameters based on ESPI measurements," in SEM Annual Conference and Exposition on Experimental and Applied Mechanics (2005) 140. G. Silva, Le Riche, R., Molimard, J., Vautrin, A., and Galerne, C., "Identification of Material Properties Using FEMU: Application to the Open Hole Tensile Test," Applied Mechanics and Materials 7 8 (2007) 141. W. J. Noh, "Mixed Mode Interfacial Fracture Toughness of Sandwich Composites at Cryogenic Temperatures," in Mechanical and Aerospace Engineering, University of Florida, Gainesville, FL (2004). PAGE 144 144 BIOGRAPHICAL SKETCH Weiqi Yin was born in Hunan Province of China in 1980. He received both his m aster s and bachelors degree in e ngineering m echanics from Tsinghua University, China in 2004 and 2002 respectively. H e also got a Master of Science in m echanical e ngineering at University of Florida in 2007. He received his Pd.D. in mechanical engineeri ng at University of Florida in the fall of 2009. Yins research interests are focused on experimental stress analysis including moir interferometry, digital image correlation, strain gage, etc., composite materials, finite element analysis and digital i mage processing 