CN101986098B - Tricolor projection-based Fourier transform three-dimensional measuring method - Google Patents

Abstract

The invention discloses a new tricolor raster projection-based Fourier transform three-dimensional measuring method, which mainly aims to accurately work out the phase distribution of raster images, and then obtain the three-dimensional topography information of objects according to the phase distribution. The method is implemented through the steps of setting the gray values of R and G components in a projected raster to follow a sinusoidal variation regulation with two different frequencies, and setting a B component as the average value of the R component or the G component, and forming a color raster image to be projected on an object to be measured; and separating the three components of an acquired deformed color raster image, then processing the separated three gray images, meanwhile, restraining background components and high-frequency noises when the Fourier transform is used to calculate relative phases, and carrying out comparison and calibration on the three gray images with the relative phase values under the two frequencies when the phases are unwrapped so as to achieve the purpose of precisely calculating absolute phases containing height information. Only one raster image is required to be projected in the whole measuring process, therefore, the method of the invention is good in real-time.

Description

Fourier transform three-dimensional measurement method based on three coloured light projections

Technical field

The invention belongs to the field of three-dimensional information reconstruct.Based on the chromatic grating projection, adopt fourier transform method to analyze the deformed grating image and find the solution relative phase, adopt dual-frequency method to carry out phase unwrapping, obtain the process of accurate absolute phase.

Background technology

Three-dimensional feature that the three-dimensional measurement technology can be described object and the three-dimensional information that obtains body surface have a wide range of applications at product detection and machining control, medical field, historical relic's protection field, aerospace field, culture field etc.In numerous three-dimensional measurement technology, the optical three-dimensional measurement technology becomes prevailing three-dimensional measurement technology with characteristics such as its non-cpntact measurement, real-time are better.

The optical three-dimensional measurement technology is to be the basis with the contemporary optics, melts optoelectronics, Computer Image Processing, and graphics, science such as signal Processing are the modern surveying technology of one.The entity that it detects optical imagery and store as information obtains the three-dimensional measurement information needed from image.With respect to other three-dimensional information measuring method, high-speed, high-precision characteristics that optical measuring technique has since the sixties in last century, have obtained extensive studies and application.Three-dimensional measurement technology based on the encode grating projection is the most representative, uses the most extensive.

Three-dimensional measurement based on optical grating projection projects to the measured object surface with raster pattern exactly, obtains the raster image of distortion by video camera, and determines the elevation information on testee relative reference plane by phase place and relation of height.When optical grating projection was to body surface, periodically the phase place of grating just received the modulation of body surface height profile, formed deformed grating, had the three-dimensional information of object in the deformed grating.The PHASE DISTRIBUTION of accurately obtaining the deformed grating image plays key effect in whole three-dimensional measurement process.

The method of obtaining PHASE DISTRIBUTION mainly contains phase-shift method, fourier transform method, wavelet analysis method.The Fourier transform mensuration is a kind of contactless method for three-dimensional measurement commonly used; It is through carrying out overall Fourier transform to an amplitude variation shape raster image; Extraction comprises the fundamental component of phase information, inverse Fourier transform, and then the relative phase value of acquisition raster image.Than traditional phase-shifting technique, only need to gather an amplitude grating image in the fourier methods measuring process, real-time is high, the real measurement that realizes dynamic object.But because the useful information that in Fourier transform, really needs to extract is a fundamental component, the frequency spectrum of background component and high fdrequency component can be overlapped with the frequency spectrum of fundamental component, causes spectral aliasing, had a strong impact on the measuring accuracy of Fourier transform three-dimensional measurement method.The measure of the inhibition spectral aliasing that adopts usually is to take the stripe pattern that two width of cloth have the π phase shift, through subtracting each other of two width of cloth images, can suppress background component effectively.But, influenced measured real-time property owing to need to gather two width of cloth images.

Because the phase value that fourier transform method obtains obtains absolute phase values between 0-2 π so need carry out phase unwrapping.But because Fourier transform has in the process that the phase place of height saltus step and relatively gentle slope object finds the solution in solution error can appear; The error of any can influence the phase unwrapping precision of next pending point when phase unwrapping; So produce continuous mistake, i.e. backguy phenomenon in the time of can causing phase unwrapping.Directly launch the scanning line method rapid speed of phase place, but robustness a little less than, very easily produce the backguy phenomenon.Dual-frequency method contains the sinusoidal grating of two kinds of frequency components through projection; Through Fourier transform filtering; Obtain two kinds of phase values under the frequency, with the phase value that obtains under the phase value correction high frequency that obtains under the low frequency, thus the improper value at the discontinuous place of phase place that obtains under the correction high frequency.But because dual-frequency method also is to find the solution initial phase with Fourier transform, so also exist two kinds of problem between the frequency component.

Summary of the invention

Technical matters:, the objective of the invention is to not influence the influence of eliminating the spectral aliasing in the Fourier transform three-dimensional measurement method under the prerequisite of measuring real-time and the precision that improves the phase unwrapping process to the spectral aliasing and the not high problem of phase unwrapping process robustness that exist in the Fourier transform three-dimensional measurement method.This method only needs projection one amplitude grating image, has improved measured real-time property, can realize kinetic measurement, and can overcome factors such as testee surface reflectivity difference and surface noise, reaches measuring accuracy preferably, has robustness preferably.

Technical scheme: a kind of measuring three-dimensional morphology method based on tricolor grating projection and Fourier transform, concrete steps are following:

Step 1: design a width of cloth colour projection raster image, and raster image is projected to the object under test surface, said colour projection raster image is provided with as follows:

f(m，n)＝a _r+b _rcos(2πf ₀m)+a _g+b _gcos(4πf ₀m)+a _b

Wherein, a _r, a _gBe R, the background light intensity of G component, b _r, b _gBe R, the reflectivity of G component, a _bBe set to R, the mean value of G two components, f ₀Be the frequency of R component sine streak, (frequency of G component sine streak is 2 times of R component sine streak for m, the n) two-dimensional coordinate of representation raster image,

Step 2: use colored CCD that body surface to be measured is taken, obtain an amplitude variation shape grating fringe image:

g(x，y)＝r _r(x，y)[a _r+b _rcos(2πf ₀x+φ ₁(x，y))]+

r _g(x，y)[a _g+b _gcos(4πf ₀x+φ ₂(x，y))]+r _b(x，y)a _b

Wherein, r _r(x, y), r _b(x, y), r _g(x, y) be object to different colours reflection of light rate, reflectivity is a constant, φ ₁(x, y), φ ₂(x y) is the PHASE DISTRIBUTION of low frequency and high frequency to be asked, (and x, the y) two-dimensional coordinate of expression deformed grating stripe pattern,

Step 3: from the deformed grating stripe pattern, separate R, G, the B component of the coloured image obtain, and from R, G component gray-scale map, remove background component, obtain removing the R component low frequency standard grayscale figure R after the background component ₁(x, y), remove the G component high frequency standard grayscale figure G after the background component ₁(x, y),

Adopt following method from R component gray-scale map, to remove background component:

Step 3.1.1: the mean value m that at first asks for R component pixel gray-scale value _rWith variance d _r, the mean value m of B component pixel gray-scale value _bWith variance d _bUse following formula to revise each grey scale pixel value of R component:

r ₁(x，y)＝r(x，y)+m _b-m _r

$R(x,y)=[r_1(x,y)-m_{r_1}]\frac{d_b}{d_r}+{\mathrm{<figure-callout\; id="1"\; label="m\; r"\; filenames="HSA00000281376100021.png"\; state="\{\{state\}\}">m</figure-callout>}}_{r_{}}$

Wherein (x y) is the R component intensity profile after proofreading and correct, r to R ₁(x y) is intermediate variable,

Be r ₁(x, mean value y),

Contain the sinusoidal phase modulation signal R component of object height information after obtaining proofreading and correct and represent the B component of background component can be expressed as as follows:

R(x，y)＝a(x，y)+b(x，y)cos[2πf ₀x+φ ₁(x，y)]；

B(x，y)＝a(x，y)；

Step 3.1.2: to R (x, y), ((x, y) (x y) deducts in two kinds of components B, is not contained the grey-scale raster images R of background and high fdrequency component from R with B then for x, y) filtering ₁(x, y):

R ₁(x，y)＝R(x，y)-B(x，y)

＝b(x，y)cos[2πf ₀x+φ ₁(x，y)]；

Adopt following method from G component gray-scale map, to remove background component:

Step 3.2.1: the mean value m that at first asks for G component pixel gray-scale value _gWith variance d _g, the mean value m of B component pixel gray-scale value _bWith variance d _bUse following formula to revise each grey scale pixel value of G component:

g ₁(x，y)＝g(x，y)+m _b-m _g

$G(x,y)=[g_1(x,y)-m_{g_1}]\frac{d_b}{d_g}+{\mathrm{<figure-callout\; id="1"\; label="m\; g"\; filenames="HSA00000281376100021.png"\; state="\{\{state\}\}">m</figure-callout>}}_{g_{}}$

Wherein (x y) is the G component intensity profile after proofreading and correct, g to G ₁(x y) is intermediate variable,

Be g ₁(x, mean value y),

Obtaining proofreading and correct the back representative contains the double frequency sinusoidal phase modulation signal G component of object height information and represents the B component of background component can be expressed as as follows:

G(x，y)＝a(x，y)+b(x，y)cos[2πf ₀x+φ ₂(x，y)]；

B(x，y)＝a(x，y)；

Step 3.2.2: to G (x, y), ((x, y) (x y) deducts in the component B, is not contained the grey-scale raster images G of background and high fdrequency component from G with B then for x, y) filtering ₁(x, y):

G ₁(x，y)＝R(x，y)-B(x，y)

＝b(x，y)cos[2πf ₀x+φ ₂(x，y)]。

Step 4: to the R component low frequency standard grayscale figure R after the removal background component ₁(x, y), remove the G component high frequency standard grayscale figure G after the background component ₁(x y) does Fourier transform, and obtain the relative phase distribution of R component and the relative phase of G component respectively and distribute,

To low frequency standard grayscale figure R ₁(x y) does Fourier transform, and it is following to obtain the process that the relative phase of R component distributes:

Step 4.1.1: with R ₁(x y) is expressed as exponential form:

$R_1(x,y)=r_r(x,y)\underset{n=-\&infin;}{\overset{n=+\&infin;}{\mathrm{\&Sigma;}}}{A}_n\mathrm{expi}[2\mathrm{\&pi;n}{f}_{0}x+n{\mathrm{\&phi;}}_1(x,y)]$

Wherein, when n=0, A _n=0,

Expression imaginary unit,

Step 4.1.2: looking y is constant, adopts 1 dimension fast fourier transform algorithm to low frequency standard grayscale figure R ₁(x, every row y) carries out Fourier transform, and conversion process is:

$R_f(f,y)={\&Integral;}_{-\&infin;}^{+\&infin;}{R}_1(x,y)\exp (-2\mathrm{\&pi;ifx})\mathrm{dx}$

The expression formula that obtains frequency domain is:

$R_f(f,y)=\underset{n=-\&infin;}{\overset{n=+\&infin;}{\mathrm{\&Sigma;}}}{Q}_n(f-n{f}_0,y)$

Wherein, Q _n(f-nf ₀, y) be A _nr _r(x, y) exp [in φ ₁(x, y)] accordingly in the expression formula of frequency domain, f representes the variable of frequency domain,

Step 4.1.3: to Q _n(f-nf ₀, y) carry out filtering and extract Q ₁(f-f ₀, y), again to Q ₁(f-f ₀, y) carry out inverse Fourier transform, obtain containing the A of phase information ₁r _r(x, y) exp [i φ ₁(x, y)], calculate A ₁r _r(x, y) exp [i φ ₁(x, y)] angle value can obtain containing the R component relative phase value φ of object height information ₁(x, y), the phase value that obtains is between 0-2 π;

To high frequency standard grayscale figure G ₁(x y) does Fourier transform, and it is following to obtain the process that the relative phase of G component distributes:

Step 4.2.1: with G ₁(x y) is expressed as exponential form:

$G_1(x,y)=r_g(x,y)\underset{n=-\&infin;}{\overset{n=+\&infin;}{\mathrm{\&Sigma;}}}{B}_n\mathrm{expi}[2\mathrm{\&pi;n}{f}_{0}x+n{\mathrm{\&phi;}}_2(x,y)]$

Wherein, when n=0, B _n=0,

Expression imaginary unit,

Step 4.2.2: looking y is constant, adopts 1 dimension fast fourier transform algorithm to low frequency standard grayscale figure R ₁(x, every row y) carries out Fourier transform, and conversion process is:

$G_f(f,y)={\&Integral;}_{-\&infin;}^{+\&infin;}{G}_1(x,y)\exp (-2\mathrm{\&pi;ifx})\mathrm{dx}$

The expression formula that obtains frequency domain is:

$G_f(f,y)=\underset{n=-\&infin;}{\overset{n=+\&infin;}{\mathrm{\&Sigma;}}}{P}_n(f-n{f}_0,y)$

Wherein, P _n(f-nf ₀, y) be B _nr _g(x, y) exp [in φ ₁(x, y)] accordingly in the expression formula of frequency domain, f representes the variable of frequency domain,

Step 4.2.3: to P _n(f-nf ₀, y) carry out filtering and extract P ₁(f-f ₀, y), again to P ₁(f-f ₀, y) carry out inverse Fourier transform, obtain containing the B of phase information ₁r _g(x, y) exp [i φ ₂(x, y)], calculate B ₁r _g(x, y) exp [i φ ₂(x, y)] angle value can obtain containing the G component relative phase value φ of object height information ₂(x, y), the phase value that obtains is between 0-2 π;

Step 5:

Step 5.1: the relative phase value φ that adopts the R component that obtains in the direct deployment step 4 of scanning Beam Method ₁(x, y) promptly through following formula point by point scanning, in PHASE DISTRIBUTION, do not have till the saltus step greater than π:

Wherein,

Be the low frequency absolute phase after launching, φ ₁(x y) is relative phase, Δ φ ₁Be the phase difference value between adjacent 2,

Step 5.2: the relative phase value φ that adopts the G component that obtains in the direct deployment step 4 of scanning Beam Method ₂(x, y) promptly through following formula point by point scanning, in PHASE DISTRIBUTION, do not have till the saltus step greater than π:

Wherein,

Be the high frequency absolute phase after launching, φ ₂(x y) is relative phase, Δ φ ₂Be the phase difference value between adjacent 2

Step 6: calculate and accurately launch phase value,

Step 6.1: Read to get the low-frequency absolute phase

and high absolute phase

Step 6.2: traversal

is calculated difference degree between the two, promptly obtains following parameter distributions:

Ceil representes the function that rounds up in the following formula,

If Q ≠ 0; Then use Q that high frequency absolute phase

is proofreaied and correct, that is:

is the accurate expansion phase value that obtains.

Step 7: read final expansion phase result

according to classical optical grating projection by launching phase result

to object height h (x; Y) conversion formula; Finally try to achieve the three-dimensional information of Measuring Object, described conversion formula is:

Wherein, l, d are the geometric parameters of measuring system, and l is the distance of projector to measurement plane, and d is the distance of CCD camera to projector,

The expression phase changing capacity, Be the expansion phase result,

Be initial phase result, ω ₀Angular frequency for projection grating.

Beneficial effect: compared with prior art; The present invention has following advantage: at first; Make full use of the RGB3 kind component of coloured image, the sinusoidal grating of two kinds of frequencies is enrolled two kinds of components of RG, the B component is set to the mean value of RG component; Greatly reduce the number of the required projection grating image of three-dimensional measurement, improved measuring speed.Because the present invention only needs projection one amplitude grating image, can realize the quick measurement to dynamic object.Secondly; Fourier transform three-dimensional measurement method than the projection of Traditional use gray level image; The present invention is owing to use the B component to be set to the mean value of RG component in projection grating; When deformed grating is handled, can when the RG component is done Fourier transform, suppress 0 frequency component, reduce the influence of spectral aliasing measuring accuracy through separating 0 frequency component that the B component obtains RG component gray level image.At last; In the phase unwrapping process; The present invention adopts and represents the RG component relative phase value of height frequency to launch method of correcting then respectively, compares traditional scanning line method, can avoid since object contain big height saltus step and noise generation separate phase error and propagation of error phenomenon.To sum up, the present invention can realize the three-dimensional information of dynamic object is measured fast and accurately, obtains accurate absolute phase distribution plan, has good real-time performance and accuracy.

Description of drawings

Fig. 1 is the process flow diagram of whole process of the present invention.

Fig. 2 is the three-dimension measuring system that the present invention adopts.

Fig. 3 is the process flow diagram that step 3 of the present invention is separated 3 kinds of components of images acquired and removal background.

Fig. 4 is that step 6 of the present invention adopts dual-frequency method correction absolute phase to obtain the process flow diagram of precise phase value.

Fig. 5 (a) is that (x, y), Fig. 5 (b)-(d) removes the gray level image G of the RGB3 kind component of background to the colour distortion raster image g that collects when adopting the inventive method to measure the motorcycle backplate ₁(x, y), R ₁(x, y), B (x, y).

Fig. 6 (a) is that (x, y) the 400th line frequency spectrum are that (x, y) the 400th line frequency spectrum (before unfiltered) (c) is low frequency standard grayscale figure R to the B component image B that represents background component (b) to R component image R before subtracting each other ₁(x, y) the 400th line frequency spectrum (d) is high frequency standard grayscale figure G ₁(x, y) the 400th line frequency spectrum.Fig. 7 is the relative phase figure under the two kinds of frequencies of height that obtain.

Fig. 8; (a) be high frequency absolute phase after launching; (b) be high frequency absolute phase

after launching

Fig. 9 is the absolute phase figure

that adopts the inventive method finally to obtain

Figure 10 is the motorcycle backplate three-dimensional point cloud design sketch that finally obtains.

Figure 11 is the geometric representation of the three-dimension measuring system that adopts of the present invention.

Embodiment

Further specify below in conjunction with the accompanying drawing specific embodiments of the invention.Under windows operating system, select for use VC++6.0, the deformed grating image that CCD collects is handled as programming tool.This instance adopts the motorcycle backplate as testee, finally obtains the more accurate whole audience absolute phase that contains the backplate three-dimensional information and distributes.

Fig. 1 is the process flow diagram of whole process of the present invention.

Fig. 2 is the three-dimension measuring system that the present invention adopts.With projector the tricolor grating image projection is arrived testee, adopt colored CCD to gather the deformed grating image, transfer in the computing machine through graph card and handle.

Fig. 3 is the process flow diagram that step 3 of the present invention is separated 3 kinds of components of images acquired and removal background.

Fig. 4 is that step 6 of the present invention adopts dual-frequency method correction absolute phase to obtain the process flow diagram of precise phase value.

The present invention utilizes three kinds of components of RGB of coloured image to solve spectral aliasing and the Phase unwrapping that exists in the Fourier transform mensuration.With two kinds of components of sinusoidal variations rule modulation RG of two kinds of frequencies, the B component is set to the mean value of the gray scale of R, G component, forms a color image frame and projects body surface.Use CCD to gather the deformed grating image that receives the modulation of testee height.Through separating 3 kinds of components in this image, B component gray level image is deducted from RG component gray level image respectively, adopt Fourier transform analysis to find the solution then and obtain two kinds of relative phase values under the frequency.Phase unwrapping is carried out in distribution to relative phase, at first adopts scanning line method to launch the high and low frequency absolute phase respectively, proofreaies and correct the high frequency absolute phase to eliminate error with the low frequency absolute phase then, obtains the higher whole audience absolute phase values of precision.Then according to the phase place of classical optical grating projection to the height conversion formula, finally obtain the three-dimensional information of Measuring Object.

The present invention is based on the Fourier transform three-dimensional measurement method of three coloured light projections, detailed step is following:

Step 1: design a width of cloth colour projection raster image, and raster image is projected to the object under test surface, said colour projection raster image is provided with as follows:

f(m，n)＝a _r+b _rcos(2πf ₀m)+a _g+b _gcos(4πf ₀m)+a _b

Wherein, a _r, a _gBe R, the background light intensity of G component, b _r, b _gBe R, the reflectivity of G component, a _bBe set to R, the mean value of G two components, f ₀Be the frequency of R component sine streak, (frequency of G component sine streak is 2 times of R component sine streak for m, the n) two-dimensional coordinate of representation raster image.

Step 2: use colored CCD that body surface to be measured is taken, obtain an amplitude variation shape grating fringe image:

g(x，y)＝r _r(x，y)[a _r+b _rcos(2πf ₀x+φ ₁(x，y))]+

r _g(x，y)[a _g+b _gcos(4πf ₀x+φ ₂(x，y))]+r _b(x，y)a _b

Wherein, r _r(x, y), r _b(x, y), r _g(x is that object is to different colours reflection of light rate, φ y) ₁(x, y), φ ₂(x y) is the PHASE DISTRIBUTION of low frequency and high frequency to be asked, (and x, the y) two-dimensional coordinate of expression deformed grating stripe pattern,

Because surround lighting is constant, so can be with reflectivity r _r(x, y), r _b(x, y), r _g(x y) is regarded as constant.

Practical distortion grating fringe image is shown in Fig. 5 (a).

Step 3: from the deformed grating stripe pattern, separate R, G, the B component of the coloured image that obtains, and from R, G component gray-scale map, remove background component, obtain low frequency standard grayscale figure R ₁(x, y), high frequency standard grayscale figure G ₁(x, y),

3 kinds of components of RGB in the deformed grating image that receives testee height modulation that CCD is collected separate; Because object is different to the reflectivity of every kind of component of RGB; So need the value of averaging with variance correction make every kind of component background component contrast identical so that do further processing.With the B component is benchmark, with the mean value of G component and R component and contrast correction to identical with the B component.

Adopt following method from R component gray-scale map, to remove background component:

Step 3.1.1: the mean value m that at first asks for R component pixel gray-scale value _rWith variance d _r, the mean value m of B component pixel gray-scale value _bWith variance d _bUse following formula to revise each grey scale pixel value of R component:

r ₁(x，y)＝r(x，y)+m _b-m _r

$R(x,y)=[r_1(x,y)-m_{r_1}]\frac{d_b}{d_r}+m_{{\mathrm{<figure-callout\; id="1"\; label="r"\; filenames="HSA00000281376100021.png"\; state="\{\{state\}\}">r</figure-callout>}}_1}$

Wherein (x y) is the R component intensity profile after proofreading and correct, r to R ₁(x y) is intermediate variable, Be r ₁(x, mean value y),

Contain the sinusoidal phase modulation signal R component of object height information after obtaining proofreading and correct and represent the B component of background component can be expressed as as follows:

R(x，y)＝a(x，y)+b(x，y)cos[2πf ₀x+φ ₁(x，y)]；

B(x，y)＝a(x，y)；

Step 3.1.2: to R (x, y), ((x, y) (x y) deducts in two kinds of components B, is not contained the grey-scale raster images R of background and high fdrequency component from R with B then for x, y) filtering ₁(x, y):

R ₁(x，y)＝R(x，y)-B(x，y)

＝b(x，y)cos[2πf ₀x+φ ₁(x，y)]；

Adopt following method from G component gray-scale map, to remove background component:

Step 3.2.1: the mean value m that at first asks for G component pixel gray-scale value _gWith variance d _g, the mean value m of B component pixel gray-scale value _bWith variance d _bUse following formula to revise each grey scale pixel value of G component:

g ₁(x，y)＝g(x，y)+m _b-m _g

$G(x,y)=[g_1(x,y)-m_{g_1}]\frac{d_b}{d_g}+m_{{\mathrm{<figure-callout\; id="1"\; label="g"\; filenames="HSA00000281376100021.png"\; state="\{\{state\}\}">g</figure-callout>}}_1}$

Wherein (x y) is the G component intensity profile after proofreading and correct, g to G ₁(x y) is intermediate variable, Be g ₁(x, mean value y),

Obtaining proofreading and correct the back representative contains the double frequency sinusoidal phase modulation signal G component of object height information and represents the B component of background component can be expressed as as follows:

G(x，y)＝a(x，y)+b(x，y)cos[2πf ₀x+φ ₂(x，y)]；

B(x，y)＝a(x，y)；

Step 3.2.2: to G (x, y), ((x, y) (x y) deducts in the component B, is not contained the grey-scale raster images G of background and high fdrequency component from G with B then for x, y) filtering ₁(x, y):

G ₁(x，y)＝R(x，y)-B(x，y)

＝b(x，y)cos[2πf ₀x+φ ₂(x，y)]°

Actual high frequency standard grayscale figure G ₁(x, y), low frequency standard grayscale figure R ₁(x, y) as Fig. 5 (b) (c) shown in.Fig. 5 (d) for the B of representative background component (x, y).Fig. 6 (a) is that (x, y) the 400th line frequency spectrum are that (x, y) the 400th line frequency spectrum (before unfiltered) (c) is low frequency standard grayscale figure R to the B component image B that represents background component (b) to R component image R before subtracting each other ₁(x, y) the 400th line frequency spectrum (d) is high frequency standard grayscale figure G ₁(x, y) the 400th line frequency spectrum.Can find out by spectrogram, through the processing of step 3, low frequency standard grayscale figure R ₁(x, y) the 400th line frequency spectrum and high frequency standard grayscale figure G ₁(x, y) background component (i.e. near frequency component 0 frequency) is suppressed fully in the 400th line frequency spectrum.

Step 4: to low frequency standard grayscale figure R ₁(x, y), high frequency standard grayscale figure G ₁(x y) does Fourier transform, and obtain the relative phase distribution of R component and the relative phase of G component respectively and distribute,

To low frequency standard grayscale figure R ₁(x y) does Fourier transform, and it is following to obtain the process that the relative phase of R component distributes:

Step 4.1.1: with R ₁(x y) is expressed as exponential form:

$R_1(x,y)=r_r(x,y)\underset{n=-\&infin;}{\overset{n=+\&infin;}{\mathrm{\&Sigma;}}}{A}_n\mathrm{expi}[2\mathrm{\&pi;n}{f}_{0}x+n{\mathrm{\&phi;}}_1(x,y)]$

Wherein, when n=0, A _n=0, Expression imaginary unit,

Step 4.1.2: looking y is constant, adopts 1 dimension fast fourier transform algorithm to low frequency standard grayscale figure R ₁(x, every row y) carries out Fourier transform, and conversion process is:

$R_f(f,y)={\&Integral;}_{-\&infin;}^{+\&infin;}{R}_1(x,y)\exp (-2\mathrm{\&pi;ifx})\mathrm{dx}$

The expression formula that obtains frequency domain is:

$R_f(f,y)=\underset{n=-\&infin;}{\overset{n=+\&infin;}{\mathrm{\&Sigma;}}}{Q}_n(f-n{f}_0,y)$

Wherein, Q _n(f-nf ₀, y) be A _nr _r(x, y) exp [in φ ₁(x, y)] accordingly in the expression formula of frequency domain, f representes the variable of frequency domain,

Step 4.1.3: to Q _n(f-nf ₀, y) carry out filtering and extract Q ₁(f-f ₀, y), again to Q ₁(f-f ₀, y) carry out inverse Fourier transform, obtain containing the A of phase information ₁r _r(x, y) exp [i φ ₁(x, y)], calculate A ₁r _r(x, y) exp [i φ ₁(x, y)] angle value can obtain containing the R component relative phase value φ of object height information ₁(x, y), the phase value that obtains is between 0-2 π;

To high frequency standard grayscale figure G ₁(x y) does Fourier transform, and it is following to obtain the process that the relative phase of G component distributes:

Step 4.2.1: with G ₁(x y) is expressed as exponential form:

$G_1(x,y)=r_g(x,y)\underset{n=-\&infin;}{\overset{n=+\&infin;}{\mathrm{\&Sigma;}}}{B}_n\mathrm{expi}[2\mathrm{\&pi;n}{f}_{0}x+n{\mathrm{\&phi;}}_2(x,y)]$

Wherein, when n=0, B _n=0,

Expression imaginary unit,

Step 4.2.2: looking y is constant, adopts 1 dimension fast fourier transform algorithm to low frequency standard grayscale figure R ₁(x, every row y) carries out Fourier transform, and conversion process is:

$G_f(f,y)={\&Integral;}_{-\&infin;}^{+\&infin;}{G}_1(x,y)\exp (-2\mathrm{\&pi;ifx})\mathrm{dx}$

The expression formula that obtains frequency domain is:

$G_f(f,y)=\underset{n=-\&infin;}{\overset{n=+\&infin;}{\mathrm{\&Sigma;}}}{P}_n(f-n{f}_0,y)$

Wherein, P _n(f-nf ₀, y) be B _nr _g(x, y) exp [in φ ₁(x, y)] accordingly in the expression formula of frequency domain, f representes the variable of frequency domain,

Step 4.2.3: to P _n(f-nf ₀, y) carry out filtering and extract P ₁(f-f ₀, y), again to P ₁(f-f ₀, y) carry out inverse Fourier transform, obtain containing the B of phase information ₁r _g(x, y) exp [i φ ₂(x, y)], calculate B ₁r _g(x, y) exp [i φ ₂(x, y)] angle value can obtain containing the G component relative phase value φ of object height information ₂(x, y), the phase value that obtains is between 0-2 π;

Fig. 7 is the R component relative phase value φ that obtains ₁(x, y), G component relative phase value φ ₂(x, relative phase figure y);

Step 5:

Step 5.1: the relative phase value φ that adopts the R component that obtains in the direct deployment step 4 of scanning Beam Method ₁(x, y) promptly through following formula point by point scanning, in PHASE DISTRIBUTION, do not have till the saltus step greater than π:

Wherein,

Be the low frequency absolute phase after launching, φ ₁(x y) is relative phase, Δ φ ₁Be the phase difference value between adjacent 2,

Step 5.2: the relative phase value φ that adopts the G component that obtains in the direct deployment step 4 of scanning Beam Method ₂(x, y) promptly through following formula point by point scanning, in PHASE DISTRIBUTION, do not have till the saltus step greater than π:

Wherein,

Be the high frequency absolute phase after launching, φ ₂(x y) is relative phase, Δ φ ₂Be the phase difference value between adjacent 2

Fig. 8; (a) be high frequency absolute phase

after launching; (b) be high frequency absolute phase

after launching

Step 6: by the absolute phase values under two kinds of frequencies that obtain in the step 5, adopt dual-frequency method to proofread and correct the high frequency absolute phase values, eliminate error, obtain final expansion phase result.

Adopt scanning line method in the step 5,, can make can make mistakes in the absolute phase values subregion under the high frequency because the shade of some sheer face can block the mistake transmission of whole cycle and big noise spot.Obtain absolute phase value under the low frequency because fringe period is enough big; The phase unwrapping mistake that shade produces is able to avoid; But because striped underfrequency; Measuring accuracy is the absolute phase values as obtaining under the high frequency scarcely, so the absolute phase values

that adopts the absolute phase values

under the low frequency to proofread and correct under the high frequency can be eliminated the zone errors in

when not reducing the phase place solving precision.

Step 6.1: Read to get the low-frequency absolute phase

and high absolute phase

Step 6.2: traversal is calculated difference degree between the two, promptly obtains following parameter distributions:

Ceil representes the function that rounds up in the following formula,

If Q ≠ 0; Then use Q that high frequency absolute phase

is proofreaied and correct, that is:

is the accurate expansion phase value that obtains.Phase diagram is as shown in Figure 9.

Step 7: read final expansion phase result

according to the phase place of classical optical grating projection to the height conversion formula, finally try to achieve the three-dimensional information of Measuring Object.

The geometric representation of measuring system is shown in figure 11; By launch phase result

to object height h (x, formula y) is following:

Wherein, l, d are the geometric parameters of measuring system, and l is the distance of projector to measurement plane, and d is the distance of CCD camera to projector,

The expression phase changing capacity,

Be the expansion phase result,

Be the initial phase result, by the decision of witness mark face, ω ₀Be the angular frequency of projection grating, can get by system calibrating.

The point cloud chart of expression motorcycle backplate three-dimensional information is shown in figure 10.

Claims (3)

1. measuring three-dimensional morphology method based on tricolor grating projection and Fourier transform is characterized in that concrete steps are following:

Step 1: design a width of cloth colour projection raster image, and raster image is projected to the object under test surface, said colour projection raster image is provided with as follows:

f(m，n)＝a _r+b _rcos(2πf ₀m)+a _g+b _gcos(4πf ₀m)+a _b

Wherein, a _r, a _gBe R, the background light intensity of G component, b _r, b _gBe R, the reflectivity of G component, a _bBe set to R, the mean value of the background light intensity of G two components, f ₀Be the frequency of R component sine streak, (frequency of G component sine streak is 2 times of R component sine streak, step 2 for m, the n) two-dimensional coordinate of representation raster image: use colored CCD that body surface to be measured is taken, obtain an amplitude variation shape grating fringe image:

g(x，y)＝r _r(x，y)[a _r+b _rcos(2πf ₀x+φ ₁(x，y))]+

r _g(x，y)[a _g+b _gcos(4πf ₀x+φ ₂(x，y))]+r _b(x，y)a _b

Wherein, r _r(x, y), r _b(x, y), r _g(x, y) be object to different colours reflection of light rate, reflectivity is a constant, φ ₁(x, y), φ ₂(x y) is the PHASE DISTRIBUTION of low frequency and high frequency to be asked, (and x, the y) two-dimensional coordinate of expression deformed grating stripe pattern,

Step 3: from the deformed grating stripe pattern, separate R, G, the B component of the coloured image obtain, and from R, G component gray-scale map, remove background component, obtain removing the R component low frequency standard grayscale figure R after the background component ₁(x, y), remove the G component high frequency standard grayscale figure G after the background component ₁(x, y),

Step 4: to R component low frequency standard grayscale figure R ₁(x, y), G component high frequency standard grayscale figure G ₁(x y) does Fourier transform, obtains the relative phase distribution of R component and the relative phase of G component respectively and distributes, to the R component low frequency standard grayscale figure R after the removal background component ₁(x y) does Fourier transform, and it is following to obtain the process that the relative phase of R component distributes:

Step 4.1.1: with R ₁(x y) is expressed as exponential form:

Wherein, when n=0, A _n=0,

Expression imaginary unit,

Step 4.1.2: looking y is constant, adopts 1 dimension fast fourier transform algorithm to low frequency standard grayscale figure R ₁(x, every row y) carries out Fourier transform, and conversion process is:

The expression formula that obtains frequency domain is:

Wherein, Q _n(f-nf ₀, y) be A _nr _r(x, y) exp [in φ ₁(x, y)] accordingly in the expression formula of frequency domain, f representes the variable of frequency domain,

Step 4.1.3: to Q _n(f-nf ₀, y) carry out filtering and extract Q ₁(f-f ₀, y), again to Q ₁(f-f ₀, y) carry out inverse Fourier transform, obtain containing the A of phase information ₁r _r(x, y) exp [i φ ₁(x, y)], calculate A ₁r _r(x, y) exp [i φ ₁(x, y)] angle value can obtain containing the R component relative phase value φ of object height information ₁(x, y), the phase value that obtains is between 0-2 π;

To the G component high frequency standard grayscale figure G after the removal background component ₁(x y) does Fourier transform, and it is following to obtain the process that the relative phase of G component distributes:

Step 4.2.1: with G ₁(x y) is expressed as exponential form:

Wherein, when n=0, B _n=0, Expression imaginary unit,

Step 4.2.2: looking y is constant, adopts 1 dimension fast fourier transform algorithm to low frequency standard grayscale figure R ₁(x, every row y) carries out Fourier transform, and conversion process is:

The expression formula that obtains frequency domain is:

Wherein, P _n(f-nf ₀, y) be B _nr _g(x, y) exp [in φ ₁(x, y)] accordingly in the expression formula of frequency domain, f representes the variable of frequency domain,

Step 4.2.3: to P _n(f-nf ₀, y) carry out filtering and extract P ₁(f-f ₀, y), again to P ₁(f-f ₀, y) carry out inverse Fourier transform, obtain containing the B of phase information ₁r _g(x, y) exp [i φ ₂(x, y)], calculate B ₁r _g(x, y) exp [i φ ₂(x, y)] angle value can obtain containing the G component relative phase value φ of object height information ₂(x, y), the phase value that obtains is between 0-2 π;

Step 5:

Step 5.1: the relative phase value φ that adopts the R component that obtains in the direct deployment step 4 of scanning Beam Method ₁(x, y) promptly through following formula point by point scanning, in PHASE DISTRIBUTION, do not have till the saltus step greater than π:

Wherein,

Be the low frequency absolute phase after launching, φ ₁(x y) is relative phase, Δ φ ₁Be the phase difference value between adjacent 2,

Step 5.2: the relative phase value φ that adopts the G component that obtains in the direct deployment step 4 of scanning Beam Method ₂(x, y) promptly through following formula point by point scanning, in PHASE DISTRIBUTION, do not have till the saltus step greater than π:

Wherein,

Be the high frequency absolute phase after launching, φ ₂(x y) is relative phase, Δ φ ₂Be the phase difference value between adjacent 2,

Step 6: adopt dual-frequency method to proofread and correct high frequency absolute phase high frequency absolute phase

by low frequency absolute phase that obtains in the step 5

and high frequency absolute phase and eliminate error, obtain final expansion phase result

Step 7: read final expansion phase result

according to classical optical grating projection by launching phase result

to object height h (x; Y) conversion formula; Finally try to achieve the three-dimensional information of Measuring Object, described conversion formula is:

Wherein, l, d are the geometric parameters of measuring system, and l is the distance of projector to measurement plane, and d is the distance of CCD camera to projector, The expression phase changing capacity,

Be the expansion phase result,

Be initial phase result, ω ₀Angular frequency for projection grating.

2. the measuring three-dimensional morphology method based on chromatic grating projection and Fourier transform according to claim 1 is characterized in that,

Adopt following method from R component gray-scale map, to remove background component:

Step 3.1.1: the mean value m that at first asks for R component pixel gray-scale value _rWith variance d _r, the mean value m of B component pixel gray-scale value _bWith variance d _bUse following formula to revise each grey scale pixel value of R component:

r ₁(x，y)＝r(x，y)+m _b-m _r

Wherein (x y) is the R component intensity profile after proofreading and correct, r to R ₁(x y) is intermediate variable,

Be r ₁(x, mean value y),

Contain the sinusoidal phase modulation signal R component of object height information after obtaining proofreading and correct and represent the B component of background component can be expressed as as follows:

R(x，y)＝a(x，y)+b(x，y)cos[2πf ₀x+φ ₁(x，y)]；

B(x，y)＝a(x，y)；

Step 3.1.2: to R (x, y), ((x, y) (x y) deducts in two kinds of components B, is not contained the R component low frequency standard grayscale figure R of background and high fdrequency component from R with B then for x, y) filtering ₁(x, y):

R ₁(x，y)＝R(x，y)-B(x，y)

＝b(x，y)cos[2πf ₀x+φ ₁(x，y)]；

Adopt following method from G component gray-scale map, to remove background component:

Step 3.2.1: the mean value m that at first asks for G component pixel gray-scale value _gWith variance d _g, the mean value m of B component pixel gray-scale value _bWith variance d _bUse following formula to revise each grey scale pixel value of G component:

g ₁(x，y)＝g(x，y)+m _b-m _g

Wherein (x y) is the G component intensity profile after proofreading and correct, g to G ₁(x y) is intermediate variable, Be g ₁(x, mean value y),

Obtaining proofreading and correct the back representative contains the double frequency sinusoidal phase modulation signal G component of object height information and represents the B component of background component can be expressed as as follows:

G(x，y)＝a(x，y)+b(x，y)cos[2πf ₀x+φ ₂(x，y)]；

B(x，y)＝a(x，y)；

Step 3.2.2: to G (x, y), ((x, y) (x y) deducts in the component B, is not contained the G component high frequency standard grayscale figure G of background and high fdrequency component from G with B then for x, y) filtering ₁(x, y):

G ₁(x，y)＝R(x，y)-B(x，y)

＝b(x，y)cos[2πf ₀x+φ ₂(x，y)]。

3. the measuring three-dimensional morphology method based on chromatic grating projection and Fourier transform according to claim 1 is characterized in that dual-frequency method is in the step 6:

Step 6.1: Read to get the low-frequency absolute phase?

and high absolute phase?

Step 6.2: traversal

is calculated difference degree between the two, promptly obtains following parameter distributions:

Ceil representes the function that rounds up in the following formula,

If Q ≠ 0; Then use Q that high frequency absolute phase

is proofreaied and correct, that is:

is the accurate expansion phase value that obtains.

Images