CN103713322A - Orientation higher-order coherent amplitude change-rate crack detection method - Google Patents

Abstract

The invention provides an orientation higher-order coherent amplitude change-rate crack detection method which belongs to the field of exploration of oil and gas fields. The method includes: (1) post-stack earthquake data processing: acquiring and processing earthquake data so as to obtain pre-stack earthquake data and performing overplay processing on the pre-stack earthquake data so as to obtain a post-stack earthquake data volume; (2) performing correlation correction on different orientation macro-bin seismic channels of a center point of a macro bin at a time point in the post-stack earthquake data volume obtained in the step (1); (3) calculating high-order amplitude change rate of the different orientation seismic channels and performing high-order amplitude change rate calculation on the seismic channels of the center point and seismic channels of adjacent points obtained after the correlation correction in the step (2) so as to obtain a high-order amplitude change rate of the seismic channels of the center point in four directions; (4) extracting the maximum orientation high-order amplitude change rate.

Description

A kind of orientation high-order coherent amplitude rate of change crack detection method

Technical field

The invention belongs to exploration of oil and gas field field, be specifically related to a kind of orientation high-order coherent amplitude rate of change crack detection method, for the prediction of oil field prospecting slit formation oil and gas reservoir.

Background technology

A kind of as in seismic properties of the rate of change of amplitude attribute, the most effectual way of having grown for research intensive district's band of seam hole, carbonate rock area type reservoir and karst reservoir.

Technology based on amplitude change rate is all to carry out on the basis of earthquake poststack seismic data, but traditional method is just carried out calculating conventionally between the seismic trace of an adjacent direction, this is normally inaccurate for having the FRACTURE PREDICTION of an orientation, because the unbalanced feature of the fracture development in partial points, with fuzzy characteristics in blocks, can not accurately be described by the amplitude change rate district of drawing a circle to approve in a lot of effectively district's bands.

Summary of the invention

The object of the invention is to solve the difficult problem existing in above-mentioned prior art, a kind of orientation high-order coherent amplitude rate of change crack detection method is provided, carry out the calculating of the amplitude attribute rate of change of high precision, multi-faceted, high-order, excavate orientative feature in poststack seismic data, adopt high-order more settlement method can than traditional coherent algorithm obtain higher signal to noise ratio (S/N ratio) and better noise print, therefore can obtain better fracture reservoir prediction effect.

The present invention is achieved by the following technical solutions:

An orientation high-order coherent amplitude rate of change crack detection method, comprises the following steps:

(1) poststack seismic data is processed: acquiring seismic data, and it is processed and obtains Prestack seismic data, described Prestack seismic data is carried out to overlap-add procedure and obtain post-stack seismic data body (the post-stack seismic data body in an earthquake work area is comprised of several ten thousand to millions of seismic traces);

(2) central point of the grand bin on a time point in post-stack seismic data body step (1) being obtained carries out the correlation-corrected of the grand bin seismic trace of different azimuth;

(3) calculate different azimuth seismic trace high-order amplitude change rate: the seismic trace of the consecutive point after the correlation-corrected that the seismic trace of this central point and step (2) are obtained is done the calculating of high-order amplitude change rate, obtains the high-order amplitude change rate of the seismic trace of this central point on four direction;

(4) maximum orientation high-order amplitude change rate extracts: get the maximal value in the high-order amplitude change rate of the seismic trace on the four direction that step (3) obtains, obtain the seismic trace of this central point at this temporal maximum orientation high-order amplitude change rate;

(5) repeatedly repeating step (2) obtains the maximum orientation high-order amplitude change rate of this central point on each time point to step (4), i.e. the maximum orientation high-order amplitude change rate of this central point on whole seismic trace;

(6) using whole earthquake work area other each point respectively as central point, repeating step (2) to (5) then, just obtains the maximum orientation high-order amplitude change rate of whole post-stack seismic data body.

Obtain maximum orientation high-order amplitude change rate and just can directly carry out FRACTURE PREDICTION, in conventional poststack seismic properties prediction crack method, amplitude change rate itself is exactly a kind of conventional and good FRACTURE PREDICTION attribute of effect, the present invention expands amplitude change rate method, and the calculated value that has obtained maximum orientation high-order amplitude change rate just can be indicated growth region and the feature in crack preferably.

In described step (2), adopt the grand bin of 5*5 or 3*3.

Described step (2) is specific as follows:

In post-stack seismic data body, get the central point of a grand bin on time point, adopt formula below the seismic trace of the seismic trace of this central point and consecutive point to be done to cross-correlation calculation one time in window T:

$r_{\mathrm{XY}}=\frac{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}(X_i-\stackrel{\‾}{X})(Y_i-\stackrel{\‾}{Y})}{\sqrt{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{(X_i-\stackrel{\‾}{X})}^2}\sqrt{\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{(Y_i-\stackrel{\‾}{Y})}^2}}$

In formula, the variable of putting centered by X, Y is the variable that other each points are done cross-correlation calculation, the value of statistical indicant that Xi, Yi are variable, i=1 to N,

be respectively the mean value of variable, N represents total number of variable point, total number at earthquake number strong point in window T while representing in the method.

Compared with prior art, the invention has the beneficial effects as follows: the inventive method is on the basis of traditional poststack seismic data, excavate as much as possible more azimuth information, it is embodied in the amplitude change rate attribute of high-order, can obtain more accurate results of fracture prediction.

Accompanying drawing explanation

Fig. 1 is the step block diagram of orientation of the present invention high-order coherent amplitude rate of change crack detection method.

Fig. 2 is grand bin central point and consecutive point seismic trace schematic diagram.

Fig. 3 is certain orientation, the oil field district high-order amplitude change rate attributed graph in the embodiment of the present invention.

Embodiment

Below in conjunction with accompanying drawing, the present invention is described in further detail:

The inventive method know-why used is as follows:

1. coherent amplitude rate of change

Conventionally amplitude change rate is defined as:

$\mathrm{VAR}(x,y,t)=\sqrt{{\left(\frac{\mathrm{dA}(x,y,t)}{\mathrm{dx}}\right)}^2+{\left(\frac{\mathrm{dA}(x,y,t)}{\mathrm{dy}}\right)}^2}---\left(1\right)$

Wherein A (x, y, t) is certain property value at place a bit in said three-dimensional body, and VAR (x, y, t) represents on this aspect along x, and the Vector Mode of the variable quantity in y direction represents the intensity of this attribute change.In actual applications, can consider that different Orthogonal Composite and difference gridding calculate and adjustment.In calculating of the present invention, can select differentiate difference gridding to comprise: front difference coefficient, rear difference coefficient, difference coefficient before and after center, quadrature and the various ways such as front and back difference coefficient of intersecting, also can consider the average of two kinds of difference coefficients in each direction, and ask the mean value of two kinds of orthogonal grid rate of change; Adjustment means have: along layer, utilize layer position average in certain hour window up and down while asking for, body utilizes several corresponding point of adjacent inclined road to slide by road while calculating, and eliminates the impact of zero point; Before calculating, for suppressing noise and making extremely have certain continuity, jump in the part of data and greatly indivedual value to suppress can to add various processing (three winks, road integral differential+filtering, level and smooth).

2. coherent amplitude rate of change

In order to make to calculate the impact that is not subject to stratum fluctuations, or eliminate the high value point of subregion amplitude and locate to have time shift to change the background influence of the non-reservoir quality that may cause, thereby in high-precision calculating the present invention introduced the relevant time shift of calculating between road by complete edge in length and breadth to the amplitude change rate of plane computations, along exhibition to three dimensions, consider the space-variant vector of amplitude, its mould is that the rate of change of amplitude is:

$\mathrm{VAR}(x,y,t,{\mathrm{\Δt}}_x,t_y^{\mathrm{\Δ}})=\sqrt{{\left(\frac{\mathrm{dA}(x,y,t,{\mathrm{\Δt}}_x)}{\mathrm{dx}}\right)}^2+{\left(\frac{\mathrm{dA}(x,y,t,{\mathrm{\Δt}}_y)}{\mathrm{dy}}\right)}^2}---\left(2\right)$

Δ tx wherein, Δ ty represents that difference coefficient calculates the relevant time shift amount of the wiggle trace of point-to-point transmission, and correlativity Gao road is just carried out to part time shift processing, represents that base value reflection truly has the impact of fluctuations; General related coefficient will not be carried out the processing in time shift, thereby considers the variation of some mixed and disorderly part amplitudes under mutually.Such improvement makes the calculating of amplitude change rate be more suitable for the crack area that stratum reflection has certain fluctuations.

3. high-order amplitude change rate formula

The approximate expression of the higher order accuracy of derivative is applied to (as Fornberg, 1987) among the digital simulation of geophysical phenomena conventionally.If take to improve the method for calculating exponent number in calculating, just can further improve the precision that amplitude change rate calculates.As an illustration of the high-order approximation of first order derivative, the present invention launches to release an expression formula that has the first order derivative of 4 rank precision by the Taylor of function u (x ± h) and u (x ± 2h) (wherein h=is CDP interval).By Taylor, launch to have:

$u(x+h)=u\left(x\right)+h\frac{d}{\mathrm{dx}}u\left(x\right)+\frac{{\mathrm{<figure-callout\; id="2"\; label="h"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">h</figure-callout>}}^2}{2!}\frac{{\mathrm{<figure-callout\; id="2"\; label="d"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">d</figure-callout>}}^2}{{\mathrm{dx}}^2}u\left(x\right)+\frac{{\mathrm{<figure-callout\; id="3"\; label="h"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">h</figure-callout>}}^3}{3!}\frac{{\mathrm{<figure-callout\; id="3"\; label="d"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">d</figure-callout>}}^3}{{\mathrm{dx}}^3}u\left(x\right)+\frac{{\mathrm{<figure-callout\; id="4"\; label="h"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">h</figure-callout>}}^4}{4!}\frac{{\mathrm{<figure-callout\; id="4"\; label="d"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">d</figure-callout>}}^4}{{\mathrm{dx}}^4}u\left(x\right)+O\left(h^5\right)$

$u(x-h)=u\left(x\right)-h\frac{d}{\mathrm{dx}}u\left(x\right)+\frac{{\mathrm{<figure-callout\; id="2"\; label="h"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">h</figure-callout>}}^2}{2!}\frac{{\mathrm{<figure-callout\; id="2"\; label="d"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">d</figure-callout>}}^2}{{\mathrm{dx}}^2}u\left(x\right)-\frac{{\mathrm{<figure-callout\; id="3"\; label="h"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">h</figure-callout>}}^3}{3!}\frac{{\mathrm{<figure-callout\; id="3"\; label="d"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">d</figure-callout>}}^3}{{\mathrm{dx}}^3}u\left(x\right)+\frac{{\mathrm{<figure-callout\; id="4"\; label="h"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">h</figure-callout>}}^4}{4!}\frac{{\mathrm{<figure-callout\; id="4"\; label="d"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">d</figure-callout>}}^4}{{\mathrm{dx}}^4}u\left(x\right)+O\left(h^5\right)---\left(3\right)$

$u(x+2h)=u\left(x\right)+2h\frac{d}{\mathrm{dx}}u\left(x\right)+\frac{{\left(2h\right)}^2}{2!}\frac{{\mathrm{<figure-callout\; id="2"\; label="d"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">d</figure-callout>}}^2}{{\mathrm{dx}}^2}u\left(x\right)+\frac{{\left(2h\right)}^3}{3!}\frac{{\mathrm{<figure-callout\; id="3"\; label="d"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">d</figure-callout>}}^3}{{\mathrm{dx}}^2}u\left(x\right)+\frac{{\left(2h\right)}^4}{4!}\frac{{\mathrm{<figure-callout\; id="4"\; label="d"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">d</figure-callout>}}^4}{{\mathrm{dx}}^4}u\left(x\right)+O({\mathrm{<figure-callout\; id="5"\; label="h"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">h</figure-callout>}}^5$

$u(x-2h)=u\left(x\right)-2h\frac{d}{\mathrm{dx}}u\left(x\right)+\frac{{\left(2h\right)}^2}{2!}\frac{{\mathrm{<figure-callout\; id="2"\; label="d"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">d</figure-callout>}}^2}{{\mathrm{dx}}^2}u\left(x\right)-\frac{{\left(2h\right)}^3}{3!}\frac{{\mathrm{<figure-callout\; id="3"\; label="d"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">d</figure-callout>}}^3}{{\mathrm{dx}}^2}u\left(x\right)+\frac{{\left(2h\right)}^4}{4!}\frac{{\mathrm{<figure-callout\; id="4"\; label="d"\; filenames="HDA00002217524800021.png"\; state="\{\{state\}\}">d</figure-callout>}}^4}{{\mathrm{dx}}^4}u\left(x\right)+O\left(h^5\right)$

Above-mentioned four expansions are associated with:

$\frac{1}{6}[u(x+2h)-u(x-2h)]-\frac{4}{3}[u(x+h)-u(x-h)]=2h\frac{d}{\mathrm{dx}}u\left(x\right)$

$\frac{d}{\mathrm{dx}}u\left(x\right)=\frac{\frac{1}{6}\left(u(x-2h)\right)-\frac{4}{3}\left(u(x-h)\right)+\frac{4}{3}\left(u(x+h)\right)-\frac{1}{6}\left(u(x+2h)\right)}{2h}---\left(4\right)$

4 rank expression formula filter operators of first order derivative are [1/6-4/3 0 4/3-1/6]=[0.167-1.33 0 1.33-0.167].

The step of the inventive method as shown in Figure 1, comprising:

(1) conventional processing of poststack seismic data

The conventional processing of poststack seismic data, the pre-service of data, mainly for the common seismic relative amplitude preserved processing work in early stage, mainly comprises editor, bandpass filtering, true amplitude recovery, static correction, velocity analysis, residual static correction, earth's surface amplitude uniformity compensation, deconvolution and normal moveout correction etc.

(2) correlation-corrected of different azimuth seismic trace (be correlated with and refer to that the signal of two outs of phase becomes the signal of same phase by correction)

Fig. 2 is that Qu Yige center seismic trace carrys out distance explanation, this grand bin can be 3. * 3 points, also can be 5. * 5 points, for the purpose of symmetry, generally do not select even number point (for example 2*2 or 4*4) to do grand bin and calculate, the road number that 7 grand bins relate to is too much, calculates too complex, because the elongation of track pitch, is also almost seldom used.

The schematic diagram that post-stack seismic data is as shown in Figure 2 arranged (all so arrange by all post-stack seismic datas, it is to be noted, what choose here is a little, the i.e. central point of a poststack seismic data, if the words of frontier point, likely on limit, do not have 3 roads can carry out grand bin expansion, common computing method can be cast out two lines on border, since the 3rd line, count), Fig. 2 mid point 13 is the central point of the grand bin on a time point of 5 point * 5, by this central point successively with point 1, point 3, point 5, point 5, point 7, point 8, point 9, point 11, point 12, point 13, point 14, point 15, point 17, point 18, point 19, point 21, point 23, the seismic trace of point 25 is done cross-correlation calculation in the time of one in window T, concrete correlation formula is:

$r_{\mathrm{XY}}=\frac{\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}(X_i-\stackrel{\&OverBar;}{X})(Y_i-\stackrel{\&OverBar;}{Y})}{\sqrt{\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{(X_i-\stackrel{\&OverBar;}{X})}^2}\sqrt{\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{(Y_i-\stackrel{\&OverBar;}{Y})}^2}}---\left(5\right)$

In formula, X, Y is respectively variable, Xi, the value of statistical indicant that Yi is variable (i=1,2,3 ... n), for the mean value of independent variable, in the present embodiment, X is the variate-value of point 13, and Y is the variate-value that other each points are done cross-correlation calculation.

Time window T in, time when related coefficient is maximum is designated as t.Centered by time t, put the correlation-corrected time of seismic trace and consecutive point seismic trace, by (5) formula, completed correlation-corrected.

(3) different azimuth seismic trace high-order amplitude change rate calculates (this step has been used the seismic trace of step 1, and in step 2 seismic trace after correlation-corrected)

The seismic trace of the consecutive point after the seismic trace of this time centre point and correlation-corrected is done to high-order amplitude change rate and calculate, concrete computing method are: at directions X, do amplitude change rate calculate (referring to know-why 1.) with point 11 to the seismic trace of point 15, seismic trace coefficient is successively 1/6 ,-4/3,0,4/3,-1/6 (referring to know-why 3.), in Y-direction, with 3,8, seismic trace on 13,18,23 participates in calculating; In miter angle direction, with 5,9, the seismic trace of 13,17,21 participates in computing; In 135 degree angular direction, with 1,7, the seismic trace of 13,19,25 participates in calculating.So just can obtain the value of the high-order amplitude change rate of the seismic trace of central point 13 on four direction.Specifically be calculated as follows:

At directions X, 11 to 15 5 points for being concerned with and calculating, for central point 13, orientation is now 0 degree (or 180 degree); In Y-direction, 3,8,13,18,23, for 13 of central points, orientation is now 90 degree (or 270 degree); Equally, for 5,9,13,17,21 mean that orientation is 45 degree (or 225 degree), and for 1,7,13,19,25 orientation are 135 degree (or 315 degree).In all orientation, all take central point 13 as 5 seismic traces calculate, 5 seismic traces are distributed and are labeled as 1,2,3,4,5. the index contrast of these 5 seismic traces corresponds to 1/6 ,-4/3,0,4/3 ,-1/6 in molecular moiety in formula (4).

$\frac{d}{\mathrm{dx}}u\left(x\right)=\frac{\frac{1}{6}\left(u(x-2h)\right)-\frac{4}{3}\left(u(x-h)\right)+\frac{4}{3}\left(u(x+h)\right)-\frac{1}{6}\left(u(x+2h)\right)}{2h}---\left(4\right)$

All use same formula to carry out calculating the seismic trace in 1 to 5 in all orientation, just can obtain the values that 0 degree, 45 degree, 90 degree and 135 of central point (13) spend and (be the high-order amplitude change rate in 4 orientation of this seismic trace (seismic trace has 1000 to 1000 data points not etc.) on this central point (13 point).

If use the grand bin of 3*3 herein, formula only goes the adjacent twice of central point accordingly, and formula is $\frac{d}{\mathrm{dx}}u\left(x\right)=\frac{-\frac{4}{3}\left(u(x-h)\right)+\frac{4}{3}\left(u(x+h)\right)}{2h}.$

(4) maximum orientation high-order amplitude change rate

By step (1), arrive step (3), can obtain the high-order amplitude change rate of central point seismic trace poststack seismic trace on four direction on a certain central point, get the maximal value in this four direction, just obtain the maximum orientation high-order amplitude change rate of central point seismic trace on this time point.

(5) repeating step (2), to (4), calculates each temporal value on this seismic trace successively, obtains the maximum orientation high-order amplitude change rate of central point on whole seismic trace;

(6) using whole earthquake work area other each point respectively as central point, repeating step (2) to (5) then, just obtains the maximum orientation high-order amplitude change rate of whole post-stack seismic data body.

It must be noted that: due to the restriction of bin and acquisition system, enough in order to guarantee the number of central point periphery bin seismic trace, conventionally do not calculate several lines on limit, earthquake work area.Simultaneously due to time window restriction, while calculating, do not calculate the above half o'clock window of seismic trace on the whole sampling time and the data in latter half o'clock window on whole seismic trace yet.

Fig. 3 is certain orientation, oil field district high-order amplitude change rate attributed graph, and in figure, high-order orientation amplitude change rate band indicates Liao Gai district solution fissure to grow favourable distributive province exactly.The unusual part that the abnormal indication of block shape does not change with orientation substantially, for corrosion Dong Feng district, be little wire, dendritic extremely for and orientation relevant portion, be generally fracture development district, these wire, dendritic abnormal population distribution, corrosion Dong Feng district neighbour, are effective development area in the inferior one-level crack with seam hole association as can be seen from Fig..

The present invention is directed to the growth problem of the slit formation reservoir facing in actual oil field prospecting, the method adopting is to carry out the high-order amplitude change rate calculating research with orientation for actual poststack seismic data, poststack central point seismic trace is launched Hong Mianyuandao collection correlation-corrected, the seismic trace of consecutive point in grand bin is carried out to the calculating of high-order amplitude change rate, and this step is generalized to whole work area the most at last, obtain the maximum orientation high-order amplitude change rate attribute forecast result of whole geological data, and then the fracture intensity in this work area is predicted

Technique scheme is one embodiment of the present invention, for those skilled in the art, the invention discloses on the basis of application process and principle, be easy to make various types of improvement or distortion, and be not limited only to the described method of the above-mentioned embodiment of the present invention, therefore previously described mode is just preferred, and does not have restrictive meaning.

Claims (3)

1. an orientation high-order coherent amplitude rate of change crack detection method, is characterized in that: said method comprising the steps of:

(1) poststack seismic data is processed: acquiring seismic data, and it is processed and obtains Prestack seismic data, described Prestack seismic data is carried out to overlap-add procedure and obtain post-stack seismic data body;

(2) central point of the grand bin on a time point in post-stack seismic data body step (1) being obtained carries out the correlation-corrected of the grand bin seismic trace of different azimuth;

(3) calculate different azimuth seismic trace high-order amplitude change rate: the seismic trace of the consecutive point after the correlation-corrected that the seismic trace of this central point and step (2) are obtained is done the calculating of high-order amplitude change rate, obtains the high-order amplitude change rate of the seismic trace of this central point on four direction;

(4) maximum orientation high-order amplitude change rate extracts: get the maximal value in the high-order amplitude change rate of the seismic trace on the four direction that step (3) obtains, obtain the seismic trace of this central point at this temporal maximum orientation high-order amplitude change rate;

(5) repeatedly repeating step (2) obtains the maximum orientation high-order amplitude change rate of this central point on each time point to step (4), i.e. the maximum orientation high-order amplitude change rate of this central point on whole seismic trace;

(6) using whole earthquake work area other each point respectively as central point, repeating step (2) to (5) then, just obtains the maximum orientation high-order amplitude change rate of whole post-stack seismic data body.

2. orientation according to claim 1 high-order coherent amplitude rate of change crack detection method, is characterized in that: the grand bin that adopts 5*5 or 3*3 in described step (2).

3. orientation according to claim 2 high-order coherent amplitude rate of change crack detection method, is characterized in that: specific as follows: described step (2) is specific as follows:

In post-stack seismic data body, get the central point of a grand bin on time point, adopt formula below the seismic trace of the seismic trace of this central point and consecutive point to be done to cross-correlation calculation one time in window T:

$r_{\mathrm{XY}}=\frac{\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}(X_i-\stackrel{\&OverBar;}{X})(Y_i-\stackrel{\&OverBar;}{Y})}{\sqrt{\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{(X_i-\stackrel{\&OverBar;}{X})}^2}\sqrt{\underset{i=1}{\overset{N}{\mathrm{\&Sigma;}}}{(Y_i-\stackrel{\&OverBar;}{Y})}^2}}$

In formula, the variable of putting centered by X, Y is the variable that other each points are done cross-correlation calculation, the value of statistical indicant that Xi, Yi are variable, i=1 to N, total number at earthquake number strong point in window T when N represents,

be respectively the mean value of variable.

Images