0=25

The time and space formulation
of azimuth moveout

Sergey Fomel and Biondo L. Biondi

sergey@sep.stanford.edu, biondo@sep.stanford.edu

ABSTRACT

Azimuth moveout (AMO) transforms 3-D prestack seismic data from one common azimuth and offset to different azimuths and offsets. AMO in the time-space domain is represented by a three-dimensional integral operator. The operator components are the summation path, the weighting function, and the aperture. To determine the summation path and the weighting function, we derive the AMO operator by cascading dip moveout (DMO) and inverse DMO for different azimuths in the time-space domain. To evaluate the aperture, we apply a geometric approach, defining AMO as the result of cascading prestack migration (inversion) and modeling. The aperture limitations provide a consistent description of AMO for small azimuth rotations (including zero) and justify the economic efficiency of the method.

Azimuth moveout (AMO) is by definition an operator that transforms common-azimuth common-offset seismic reflection data to different azimuths and offsets. A constructive approach to AMO was proposed by Biondi and Chemingui . According to this approach, an AMO operator is built by cascading the dip moveout (DMO) operator that transforms the input common-azimuth data to zero offset, and the inverse DMO that transforms the zero-offset data to a new offset and azimuth. Evaluating the cascade of the frequency-domain DMO and inverse DMO operators by means of the stationary phase technique produces the integral (Kirchhoff-type) 3-D AMO operator in the time-space domain.

The first part of this paper applies an analogous idea to construct the AMO operator from the time-space domain DMO and achieves the same result in a simpler way. Cascading DMO and inverse DMO allows us to evaluate the AMO operator's summation path and the corresponding weighting function. However, it is not sufficient for evaluating the third major component of the integral operator, that is, its aperture (range of integration). To solve this problem, we apply an alternative approach, that defines AMO as a cascade of 3-D migration (inversion) for particular common-azimuth and common-offset data and 3-D modeling for a different azimuth and offset. This definition resembles the viewpoint on DMO developed by Deregowski and Rocca . As with the DMO case, the migration and modeling approach reveals the physics of the AMO aperture and limits its boundaries. It is the aperture limitation that allows us to overcome the paradoxical inconsistency between 2-D and 3-D AMO operators discussed by Biondi and Chemingui . If the aperture is chosen properly, the AMO operator converges to the 2-D offset continuation limit as the azimuth rotation approaches zero. This remarkable fact supports the proof of economical efficiency of AMO in comparison with the prestack migration operator, which is known to have an unlimited aperture.

CASCADING DMO AND INVERSE DMO IN TIME-SPACE DOMAIN In this section, we present a new version of the AMO derivation. Since the entire derivation is performed in the time-space domain, it is more straightforward than the stationary phase technique developed for the same purpose by Biondi and Chemingui .

Let $P_1\left({\bf x_1},t_1;{\bf h_1}\right)$ be the input of an AMO operator (common-azimuth and common-offset seismic reflection data after normal moveout correction) and $P_2\left({\bf x_2},t_2;{\bf h_2}\right)$ be the output. Here ${\bf x_i}\; (i=1,2)$ are midpoint locations on the surface: ${\bf x_i}=\left\{x_i,y_i\right\}$, and ${\bf h_i}$ are half-offset vectors. The 3-D AMO operator has the following general form:

$$P_2\left({\bf x_2},t_2;{\bf h_2}\right) = {\bf D}_{t_2} \int... ..., h_2}\right);\, {\bf h_1}\right)\,d{\bf x_1}\;, \EQNLABEL{AMO}$$ (1)

where ${\bf D}$ is the differentiation operator $\left({\bf D}_t\equiv{d \over dt}\right)$, $t_2\,\theta_{12}$ is the summation path, and w₁₂ is the weighting function. In this section we will evaluate $\theta_{12}$ and w₁₂ using the cascade of integral 3-D DMO and inverse DMO operators in the time-space domain. The idea of this derivation originated in Biondi and Chemingui's paper , where it was applied with the frequency-domain DMO and inverse DMO operators. In the next section, we apply a new geometric approach to evaluate the AMO aperture (range of integration in AMO).

To derive AMO in the time-space domain, an integral (Kirchoff-type) DMO operator of the form

$$P_0\left({\bf x_0},t_0;{\bf 0}\right)={\bf D}_{-t_0}^{1/2}\,... ..._1;x_0,h_1}\right); {\bf h_1}\right)\,d\hat{x}_1 \EQNLABEL{DMO}$$ (2)

is cascaded with an inverse DMO of the form

$$P_2\left({\bf x_2},t_2;{\bf h_2}\right)={\bf D}_{t_2}^{1/2}\... ...x_2, h_2}\right); {\bf 0}\right)\,d\hat{x}_0\;, \EQNLABEL{IDMO}$$ (3)

where ${\bf D}_t^{1/2}$ stands for the operator of half-order differentiation (equivalent to $(i \omega)^{1/2}$ multiplication in Fourier domain), $t_0\, \theta_{10}$ and $t_2\,\theta_{02}$ are the summation paths of the DMO and inverse DMO operators ():

$$\begin{eqnarray} \theta_{10}({\bf x_1;x_0, h_1})= \left(1-{\bf \left(x_1-x_0\rig... ...ft(x_0-x_2\right)^2 \over h_2^2}\right)^{1/2}\;, \EQNLABEL{IDMOtt}\end{eqnarray}$$ (4) (5)

w₁₀ and w₀₂ are the corresponding weighting functions (amplitudes of impulse responses), $\hat{x}_1$ is the component of ${\bf x_1}$ along the ${\bf h_1}$ azimuth, and $\hat{x}_0$ is the component of ${\bf x_0}$along the ${\bf h_2}$ azimuth. Integral operators DMO and IDMO correspond to the high-frequency asymptotic (the geometrical seismic) description of the wave field. As shown by Stovas and Fomel , operator IDMO has an asymptotically equivalent form

$$P_2\left({\bf x_2},t_2;{\bf h_2}\right)=\int \tilde{w}_{02}... ..._2, h_2}\right); {\bf 0}\right)\,d\hat{x}_0\;, \EQNLABEL{tIDMO}$$ (6)

where $\tilde{w}_{02}=w_{02}\,\sqrt{\theta_{02}}$.

Both DMO and inverse DMO operate on 3-D seismic data as 2-D operators, since their apertures are defined on a line. This implies that for a given input midpoint ${\bf x_1}$, the corresponding location of ${\bf x_0}$ must belong to the line going through ${\bf x_1}$, with the azimuth defined by the input offset ${\bf h_1}$. Similarly, ${\bf x_0}$ must be on the line going through ${\bf x_2}$ with the azimuth of ${\bf h_2}$ (Figure amox12). These theoretical facts lead us to the following conclusion:

For a given pair of input and output midpoints ${\bf x_1}$ and ${\bf x_2}$of the AMO operator, the corresponding midpoint ${\bf x_0}$ on the intermediate zero-offset gather is determined by the intersection of two lines drawn through ${\bf x_1}$ and ${\bf x_2}$in the offset directions.

Applying the geometric connection among the three midpoints, we can find the cascade of the DMO and inverse DMO operators in one step. For this purpose, it is convenient to choose an orthogonal coordinate system $\{x,y\}$ on the surface in such a way that the direction of the x axis corresponds to the input azimuth (Figure amox12). In this case the connection between the three midpoints is given by

$$\begin{eqnarray} y_0=y_1\;;\;x_0=x_2-\left(y_2-y_1\right) \cot{\varphi}\;, \EQNL... ...at{x}_1=dx_1\;;\;d\hat{x}_0=dy_1 \csc{\varphi}\;. \EQNLABEL{jacob}\end{eqnarray}$$ (7) (8)

 

amox12
Figure 1 Geometric relationships between input and output midpoint locations in AMO.

Substituting DMO into tIDMO and taking into account jacob produces the 3-D integral AMO operator AMO, where

$$\begin{eqnarray} \theta_{12}\left({\bf x_1;x_2, h_2}\right) & = & \theta_{02}\le... ...\left(y_2-y_1\right)\,\cos{\varphi}\right)^2}}\;, \EQNLABEL{AMOtt}\end{eqnarray}$$ (9)

$$= w_{02}\left({\bf x_0;x_2, h_2},t_2\right)\, w_{10}\left({\... ...theta_{02}\left({\bf x_0;x_2, h_2}\right)}}\;, \EQNLABEL{AMOwf}$$ (10)

$d{\bf x_1}=dx_1\,dy_1$. Equation $t_1=t_2\,\theta_{12}\left({\bf x_1;x_2, h_2}\right)$ is the same as equation (4) in () except for a different notation. The weighting function of the derived AMO operator $\left(w_{12}\right)$ depends on the weighting functions of DMO and inverse DMO that are involved in the construction. In Appendix A, we apply equation AMOwf to two popular versions of the DMO weighting functions that correspond to Hale's and Zhang's DMO operators.

Deriving formula AMOtt, we have to assume that the input and output offset azimuths are different ($\varphi \neq 0$). In the case of equal azimuths, AMO reduces to 2-D offset continuation (OC). The location of ${\bf x_0}$ in this case is not constrained by the input and output midpoints and can take different values on the line. Therefore the superposition of DMO and inverse DMO for offset continuation is a convolution on that line. To find the summation path of the OC operator, we should consider the envelope of the family of traveltime curves (where x₀ is the parameter of a curve in the family):

$$t_1=t_2\,\theta_{12}\left(x_1;x_2, h_2\right) = t_2\,{\left\... ...2} \over {h_1^2-\left(x_1-x_0\right)^2}}\;. \EQNLABEL{OCfamily}$$ (11)

Solving the envelope condition

$${\partial \theta_{12} \over \partial x_0}=0 \EQNLABEL{OCenvelop}$$ (12)

with respect to x₀ produces

$$x_0={{\left(\Delta\,x\right)^2+h_2^2-h_1^2+ \mbox{sign}\left... ...\,h_2^2}} \over {2\,\left(\Delta\,x\right)}}\;, \EQNLABEL{OCx0}$$ (13)

where $\Delta\,x=x_1-x_2$.Substituting OCx0 into OCfamily, we get the explicit expression of the OC summation path:

$$\begin{eqnarray} t_1 & = & \nonumber \\ & & {t_2 \over \left\vert h_2\right\vert... ...,\sqrt{2 \over {U+V}} \;\mbox{for $h_2 < h_1$}\;, \EQNLABEL{OCtt}\end{eqnarray}$$ (14)

where $U=h_1^2+h_2^2-\left(\Delta\,x\right)^2$, and $V=\sqrt{U^2-4\,h_1^2h_2^2}$.Equation OCtt corresponds to formula (6) in () (with a typo corrected). The same expression was obtained in a different way by Stovas and Fomel . The apparent difference between the 2-D and 3-D solutions introduces the problem of finding a consistent description valid for both cases. Such a description is especially important for practical applications dealing with small angles of azimuth rotation, e.g. cable feather correction in marine seismics. The next section develops a way of solving this problem, which refers to the kinematic theory of AMO and follows the ideas that Deregowski and Rocca applied to DMO-type operators.

AMO APERTURE: CASCADING MIGRATION AND MODELING The impulse response of the AMO operators corresponds to a spike on the initial constant-offset constant-azimuth gather. Such a spike can physically occur in the case of a focusing ellipsoidal reflector whose focuses are coincident with the initial source and receiver locations (the impulse response of prestack common-offset migration). Therefore, the impulse response of AMO corresponds kinematically to a reflection from this ellipsoid. These considerations allow us to define AMO as the cascade of the 3-D common-offset common-azimuth migration and the 3-D modeling for a different azimuth and offset. An analogous point of view was developed for the 2-D case by Deregowski and Rocca .

Let's consider the general symmetric ellipsoid equation

$$z(x,y)=\sqrt{R^2-\beta\,\left(x-x_1\right)^2-\left(y-y_1\right)^2}\;, \EQNLABEL{ellips}$$ (15)

where z stands for the depth coordinate, R is the small semi-axis of the ellipsoid, and $\beta$ is a nondimensional parameter describing the stretching of the ellipse $(\beta < 1)$. Deregowski and Rocca derived the following connections between the geometric properties of the reflector and the coordinates of the corresponding spike in the data:

$$R={{v\,t_1}\over 2}\;;\; \beta={t_1^2 \over t_1^2+{{4\,{\bf h_1}^2}\over v^2}}\;, \EQNLABEL{rbeta}$$ (16)

where v is the wave velocity. The center of the ellipsoid is at the initial midpoint ${\bf x_1}$.

This section addresses the kinematic problem of reflection from the ellipsoid defined by ellips. In particular, we are looking for the answer to the following question: For a given elliptic reflector defined by the input midpoint, offset, and time coordinates, what points on the surface can form a source-receiver pair valid for a reflection? If a point in the output midpoint-offset space cannot be related to a reflection pattern, we should exclude it from the AMO impulse response defined in AMO.

Fermat's principle provides a general method of solving the kinematic reflection problems. Consider a formal expression for the two-point reflection traveltime

$$t={\sqrt{ (\mbox{\boldmath{$s - \xi$}})^2+z^2(\xi_x,\xi_y)} ... ... - \xi$}})^2+z^2(\xi_x,\xi_y)} \over v}\;, \EQNLABEL{reflecttt}$$ (17)

where $\xi=$$\left\{\xi_x,\xi_y\right\}$ is the vertical projection of the reflection point to the surface, ${\bf s}=\left\{s_x,s_y\right\}={\bf x_2-h_2}$is the source location, and ${\bf r}=\left\{r_x,r_y\right\}={\bf x_2+h_2}$ is the receiver location. According to Fermat's principle, the reflection ray path between two fixed points must correspond to the extremum value of the traveltime. Hence, in the vicinity of a reflected ray,

$${\partial t \over \partial \xi_x}=0\;;\; {\partial t \over \partial \xi_y}=0\;. \EQNLABEL{fermat}$$ (18)

Solving the system of equations fermat for $\xi_x$ and $\xi_y$allows us to find the reflection ray path for a given source-receiver pair on the surface. The solution is derived in Appendix B to be

$$\xi_x={{x_0-\beta\,x_1}\over{1-\beta}}\;, \EQNLABEL{x02xi}$$ (19)

$$\xi_y=y_1+\left(x_0-\xi_x\right)\,\cot{\varphi}- {{\left(y_2... ..._2^2}\sin^2{\varphi}-\left(y_2-y_1\right)^2}}\;, \EQNLABEL{xiy}$$ (20)

where x₀ has the same meaning as in the preceding section and is defined by x012.

Since the reflection point is contained inside the ellipsoid, its projection obeys the evident inequality

$$\left(\xi_y-y_1\right)^2\leq R^2-\beta\,\left(\xi_x-x_1\right)^2\;. \EQNLABEL{xiyleq}$$ (21)

It is inequality xiyleq that defines the aperture of the AMO operator.

 

amoapp
Figure 2 The AMO impulse response traveltime. Parameters: $\left\vert{\bf h1}\right\vert=1000$ m, $\left\vert{\bf h2}\right\vert=750$ m, t₁=1 sec. The top plots illustrate the case of an unrealistically low velocity (v=10 m/s); on the bottom, v=2000 m/s. On the left side the azimuth rotation $\varphi=30^{\circ}$; on the right, $\varphi=3^{\circ}$

The AMO operator's contours for different azimuth rotation angles are shown in Figure amoapp. Comparing the results for the case of an unrealistically low velocity (the top two plots in Figure amoapp) and the case of a realistic velocity (the bottom two plots) clearly demonstrates the gain in the reduction of the aperture size achieved by the aperture limitation. The gain is especially spectacular for small azimuths. When the azimuth rotation approaches zero, the area of the 3-D aperture monotonously shrinks to a line, and the limit of the traveltime of the AMO impulse response (the inverse of AMOtt) approaches the offset continuation operator OCtt (Figures amocom). This means that taking into account the aperture limitations of AMO provides a consistent description valid for small azimuth rotations including zero (the offset continuation case). Obviously, the cost of an integral operator is proportional to its size. The size of the offset continuation operator cannot extend the difference between the offsets ${\bf \left\vert\left\vert h_1\right\vert-\left\vert h_2\right\vert\right\vert}$. If we applied DMO and inverse DMO explicitly, the total size of the two operators would be about ${\bf\left\vert h_1\right\vert+\left\vert h_2\right\vert}$, which is substantially greater. This fact proves that in the case of small azimuth rotations the AMO price is less than those of not only 3-D prestack migration, but also 3-D DMO and inverse DMO combined (). Figure amoavs shows the saddle shape of the AMO operator impulse response in a 3-D AVS display.

 

amocom Figure 3 Traveltime curves of the impulse responses. The dashed lines indicate the AMO impulse response with an azimuth rotation of 3 degrees (projection on the x plane); the solid lines, the 2-D offset continuation impulse response.

 

amoavs
Figure 4 AMO impulse response traveltime in three dimensions (the AVS display). Parameters: $\left\vert{\bf h1}\right\vert=1000$ m, $\left\vert{\bf h2}\right\vert=750$ m, t₁=1 sec, v=2000 m/s, $\varphi=30^{\circ}$.

We have applied two different theoretical approaches to AMO to find a complete definition of the integral operator AMO. Biondi and Chemingui proposed cascading the DMO and inverse DMO operators to define AMO in the frequency domain. The same approach is repeated here in a simpler way by transferring the analysis to the natural time-space domain. A new contribution to the evaluation of the AMO operator follows from applying a different approach, which extends the geometric theory of DMO () to the AMO case. Cascading prestack migration and modeling allows us to evaluate the AMO operator aperture. The compactness of the AMO aperture indicates that the integral operator can be performed at a low cost and therefore promises economic benefits for its practical implementation. [SEP,EAEG,paper,GEOTLE]

A AMO AMPLITUDE The weighting function of the AMO operator can be determined from cascading the DMO and inverse DMO operators by means of equation AMOwf. In the case of Hale's DMO () and its adjoint (),

$$w_{10}\left({\bf x_1;x_0, h_1},t_0\right) = \sqrt{t_0 \over... ...rt} \over {h_2^2-\left(x_1-x_0\right)^2}}}\;, \EQNLABEL{HDMOwf}$$ (22)

$$w_{02}\left({\bf x_0;x_2, h_2},t_2\right) = \sqrt{t_2 \over... ...t} \over {h_2^2-\left(x_0-x_2\right)^2}}}\;. \EQNLABEL{HIDMOwf}$$ (23)

As follows from HDMOwf,HIDMOwf, and AMOwf,

$$\begin{eqnarray} \times\, {{{\bf \left\vert h_1\right\vert\,\left\vert h_2\right... ...{\varphi}^2-\left(y_2-y_1\right)^2\right)}}\;. \EQNLABEL{HAMOwfdf}\end{eqnarray}$$ (24)

In the case of the so-called true-amplitude DMO () and its asymptotic inverse,

$$w_{10}\left({\bf x_1;x_0, h_1},t_0\right) = \sqrt{t_0 \over... ...,\left(h_1^2-\left(x_1-x_0\right)^2\right)}\;, \EQNLABEL{DMOwf}$$ (25)

$$w_{02}\left({\bf x_0;x_2, h_2},t_2\right) = \sqrt{t_2 \over... ...rt} \over {h_2^2-\left(x_0-x_2\right)^2}}}\;. \EQNLABEL{IDMOwf}$$ (26)

Inserting DMOwf and IDMOwf into AMOwf yields

$$\begin{eqnarray} \times\, {{{\bf h_1^2}\,\sin{\varphi}^2+ \left(\left(x_2-x_1\ri... ...n{\varphi}^2-\left(y_2-y_1\right)^2\right)}}\;. \EQNLABEL{AMOwfdf}\end{eqnarray}$$ (27)

B DERIVING THE AMO APERTURE

 

amosym Figure 5 Reflection from the ellipsoid of a prestack migration impulse response (a scheme). Top: Map view. Bottom: Section of the ellipsoid with the plane drawn through the central line and the reflection point.

This appendix describes the derivation of the main formulas for the aperture evaluation that follow from the Fermat principle fermat. In order to avoid the algebraic complications of fermat, we simplify the problem by taking into account the cylindrical symmetry of the ellipsoidal reflector ellips.

Consider a plane drawn through the reflection point and the central line of the ellipsoid (the axis of the cylindrical symmetry). This plane has to contain the central (normally reflected) ray from the reflector. This conclusion follows from the fact that all the normal reflections emerge at the central line because of the cylindrical symmetry, as shown in Figure amosym. The intersection of the 3-D reflector and the plane is the 2-D ellipse

$$\hat{z}(x)=\sqrt{R^2-\beta\,\left(x-x_1\right)^2}\;. \EQNLABEL{ellips2D}$$ (28)

The connection between the emergence point of the normal ray x₀ and the x coordinate of the reflection point $\xi_x$ can be derived from the relationship evident in Figure amosym, as follows:

$$x_0=\xi_x-\hat{z}\left(\xi_x\right)\,\tan{\alpha}= \xi_x+\ha... ...t(\xi_x\right)=\xi_x\,(1-\beta)+ \beta\,x_1\,. \EQNLABEL{xi2x0}$$ (29)

Equation xi2x0 allows us to evaluate $\xi_x$ in terms of x₀ and get x02xi. The emergence point of the normal ray x₀ corresponds to the midpoint on an imaginary zero-offset section ( with a coincident source and receiver). Therefore, the location of this point is determined for given input and output midpoints in accordance with expression x012.

Obviously, the reflection point has to be inside the ellipse ellips2D. Therefore, its projection obeys the inequality

$$\left\vert\xi_x-x_1\right\vert \leq {R\over \sqrt{\beta}}\;. \EQNLABEL{xileq}$$ (30)

As follows from xileq, xi2x0, and rbeta,

$$\left\vert x_0-x_1\right\vert \leq {R\,(1-\beta)\over \sqrt{... ...r \sqrt{{{v^2\,t_1^2}\over 2}+{\bf h_1^2}}}\;. \EQNLABEL{x0leq}$$ (31)

Inequality x0leq is the known aperture limitation of the DMO operator DMO found by Deregowski and Rocca . The equality in x0leq is achieved when the reflection point is on the surface, where the reflector dip increases to 90 degrees.

Now the only unknown left in our problem is the y-coordinate of the reflection point $\xi_y$. To find this unknown, we substitute x02xi into reflecttt, choosing the convenient parameterization

$${\bf s = x_0+ h_s\;;\; r = x_0+ h_r}\;, \EQNLABEL{x02sr}$$ (32)

where ${\bf h_r-h_s}=2\,{\bf h_2}$, and ${\bf h_r+h_s}=2\ {\bf \left(x_1-x_0\right)}$ (Figure amosym). The two-point traveltime function in reflecttt transforms to the form

$$\begin{eqnarray} t & = & {\sqrt{\left(x_0-\xi_x\right)^2-\beta\,\left(x_1-\xi_x\... ...box{\boldmath{$x_0-\xi$}}\right)} \over v}\;. \EQNLABEL{reflectrt}\end{eqnarray}$$ (33)

Applying the second equation from fermat, we get a simple linear equation for $\xi_y$, which has the explicit solution xiy. From x02xi and xiy one can find the reflection point location for given midpoint and offset. To find the limits of possible output midpoint locations, we constrain the reflection point to be inside the ellipsoid ellips similarly to the way we did in two dimensions when deriving x0leq. First, let's consider the case of y₂=y₁ (the output midpoint ${\bf x_2}$ is on the line drawn through ${\bf x_1}$ in the direction of the input azimuth). In this case, combining expression xiy and inequality xiyleq produces

$$\left\vert x_0-x_1\right\vert \leq {R\,(1-\beta)\over \sqrt{\beta+\beta^2\,\cot^2{\varphi}}}\;. \EQNLABEL{x0leq3D}$$ (34)

For any azimuth rotation angle $\varphi$ less than 90 degrees, the limitation x0leq3D is smaller than that of the DMO operator x0leq. The difference increases with the decrease of the azimuth rotation, since the AMO aperture section on the line y₂=y₁ monotonously shrinks to a point x₂=x₀=x₁ when $\varphi$ approaches zero. To extend this conclusion to the whole 3-D aperture, we can find the contour of the aperture by putting the reflection point at the edge of the ellipsoid ellips, as follows:

$$\left(\xi_y-y_1\right)^2= R^2-\beta\,\left(\xi_x-x_1\right)^2 \EQNLABEL{edge}$$ (35)

and solving xiy for y₂. The aperture contour can then be defined by the system of parametric expressions

$$y_2\left(\xi_x\right) = y_1+d\left(\xi_x\right)\,\sin{\varphi}\;, \EQNLABEL{ycon}$$ (36)

$$x_2\left(\xi_x\right) = \xi_x\,(1-\beta)+\beta\,x_1+ d\left(\xi_x\right)\,\cos{\varphi}\;, \EQNLABEL{xcon}$$ (37)

where

$$d\left(\xi_x\right)= {{d_y^2+d_x^2- \sqrt{\left(d_y^2+d_x^2... ...y\,\sin{\varphi}+d_x\,\cos{\varphi} \right)}}\;, \EQNLABEL{dxi}$$ (38)

$d_x\left(\xi_x\right)=\xi_x-x_0=\beta\,\left(\xi_x-x_1\right)$,and $d_y\left(\xi_x\right)=\xi_y-y_1$ is defined by edge.