Methods for inverting electromagnetic logging measurements

A two-step inversion method for computing multi-layer subterranean formation properties includes processing gain compensated electromagnetic measurement quantities using a first inversion to compute a corresponding set of borehole corrected gain compensated measurement quantities. The first inversion includes a mathematical model of the tool and the borehole in a uniform, anisotropic formation. The set of borehole corrected gain compensated measurement quantities are then processed using a second inversion to compute multi-layer anisotropic formation properties. The second inversion includes a 1D inversion employing a point dipole model and a multi-layer formation model.

FIELD OF THE INVENTION

Disclosed embodiments relate generally to downhole electromagnetic logging methods and more particularly to methods for inverting electromagnetic logging measurements using a layered formation model.

BACKGROUND INFORMATION

The use of electromagnetic measurements in prior art downhole applications, such as logging while drilling (LWD), is well known. Such techniques may be utilized to determine a subterranean formation resistivity, which, along with formation porosity measurements, is often used to indicate the presence of hydrocarbons in the formation.

In order to accurately measure formation properties, such as the vertical and horizontal formation resistivities, the electromagnetic logging measurements may need to be corrected for borehole effects. Such effects are related to the resistivity of the drilling fluid, the wellbore diameter, and the position of the logging tool in the wellbore. However, the position of the logging tool can vary rapidly during data acquisition owing to drill string vibration while drilling. The tool vibration dynamics are commonly rapid enough that the position of the tool in the wellbore varies significantly within the data acquisition time window. With such a data record, inverting for accurate borehole and formation properties is difficult (or even impossible).

Subterranean formations commonly include multiple layers (e.g., sedimentary layers) each having distinct electrical properties. Limited by computer processing speed, it is currently impractical (or impossible) to use a model including structural details of the logging tool and the borehole as well as including multiple formation layers having different electrical properties. There is a need for an improved inversion method that accounts for logging tool vibration and multi-layer subterranean formations.

SUMMARY

A method for computing multi-layer subterranean formation properties is disclosed. An electromagnetic logging while drilling tool acquires electromagnetic voltage measurements while rotating in a subterranean borehole. The voltage measurements are processed to compute a plurality of gain compensated measurement quantities which are further processed using a first inversion to compute a corresponding set of borehole corrected gain compensated measurement quantities. The first inversion includes a mathematical model of the tool and the borehole in a uniform, anisotropic formation. The set of borehole corrected gain compensated measurement quantities are then processed using a second inversion to compute multi-layer anisotropic formation properties. The second inversion includes a 1D inversion employing a point dipole model and a multi-layer formation model.

In one example embodiment, the borehole corrected gain compensated measurement quantities are computed via processing the electromagnetic voltage measurements using a first inversion to compute borehole and formation properties. The borehole and formation properties are further processed using the mathematical model to compute a first set of modeled gain compensated measurement quantities. The formation properties are further processed using a point dipole model to compute a second set of modeled gain compensated measurement quantities. The borehole effects may then be computed by processing differences between the first and second sets of modeled gain compensated measurement quantities. The borehole corrected gain compensated measurement quantities may be computed by subtracting the borehole effects from the gain compensated measurement quantities computed from the acquired voltage measurements.

The disclosed embodiments may provide various technical advantages. For example, the disclosed two-step inversion methodology enables multi-layer formation properties to be computed using the limitations of current computer processors (those of ordinary skill will readily appreciate that solving a single-step inversion including a multi-layer anisotropic formation as well physical details regarding the tool and borehole would be impractical or impossible due to the processing limitations of modern computers).

The disclosed embodiments may be further advantageous in that the borehole effects are removed using a point dipole model that essentially matches the 1D inversion model used compute the multi-layer formation properties. Such “matching” of the point dipole models tends to improve accuracy and reliability of the computed multi-layer formation properties. Certain disclosed embodiments further make use of gain compensated measurement quantities thereby eliminating errors due to drift of the gains in downhole environment and also eliminating expensive (and cumbersome) procedures to account for the gains.

DETAILED DESCRIPTION

FIG. 1depicts an example drilling rig10suitable for employing various method embodiments disclosed herein. A semisubmersible drilling platform12is positioned over an oil or gas formation (not shown) disposed below the sea floor16. A subsea conduit18extends from deck20of platform12to a wellhead installation22. The platform may include a derrick and a hoisting apparatus for raising and lowering a drill string30, which, as shown, extends into borehole40and includes a drill bit32deployed at the lower end of a bottom hole assembly (BHA) that further includes an electromagnetic measurement tool50configured to make directional electromagnetic logging measurements. As described in more detail below the electromagnetic measurement tool50may include multi-axial antennas deployed on a logging while drilling tool body.

It will be understood that the deployment illustrated onFIG. 1is merely an example. Drill string30may include substantially any suitable downhole tool components, for example, including a steering tool such as a rotary steerable tool, a downhole telemetry system, and one or more MWD or LWD tools including various sensors for sensing downhole characteristics of the borehole and the surrounding formation. The disclosed embodiments are by no means limited to any particular drill string configuration.

It will be further understood that the disclosed embodiments are not limited to use with a semisubmersible platform12as illustrated onFIG. 1. The disclosed embodiments are equally well suited for use with either onshore or offshore subterranean operations.

FIG. 2Adepicts one example of an electromagnetic measurement tool50. In the depicted embodiment measurement tool50includes first and second axially spaced transmitters52and54and first and second axially spaced receivers56and58deployed on a logging while drilling tool body51, with the receivers56and58being deployed axially between the transmitters52and54. To obtain directional measurements, each of the transmitters52and54and receivers56and58generally includes at least one transverse antenna and may further include an axial antenna. For example, the transmitters and receivers may include a bi-axial antenna arrangement including an axial antenna and a transverse (cross-axial) antenna. In another embodiment, the transmitters and receivers may include a tri-axial antenna arrangement including an axial antenna and first and second transverse antennas that are orthogonal to one another. As is known to those of ordinary skill in the art, an axial antenna is one whose moment is substantially parallel with the longitudinal axis of the tool. Axial antennas are commonly wound about the circumference of the logging tool such that the plane of the antenna is substantially orthogonal to the tool axis. A transverse antenna is one whose moment is substantially perpendicular to the longitudinal axis of the tool. A transverse antenna may include, for example, a saddle coil (e.g., as disclosed in U.S. Patent Publications 2011/0074427 and 2011/0238312 each of which is incorporated by reference herein).

While not depicted onFIGS. 2A-2C, it will be understood that one or more of the transmitters52and54and the receivers56and58may include a tilted antenna. Tilted antennas are commonly used to make directional resistivity measurements. As is known to those of ordinary skill in the art, a tilted antenna is one whose moment is angularly offset (tilted) with respect to the tool axis and is neither parallel with nor orthogonal to the tool axis.

FIG. 2Bdepicts the moments (magnetic dipoles) of one embodiment of measurement tool50in which the transmitters52,54and receivers56,58each include a tri-axial antenna arrangement. Each of the transmitters52,54includes an axial transmitting antenna T1zand T2zand first and second transverse transmitting antennas T1x, T1yand T2x, T2y. Likewise, each of the receivers56,58includes an axial receiving antenna R1zand R2zand first and second transverse receiving antennas R1x, R1yand R2x, R2y. It will be understood that the disclosed embodiments are not limited to a tri-axial antenna configuration such as that depicted onFIG. 2B.

FIG. 2Cdepicts an alternative electromagnetic measurement tool embodiment50′ in which the first and second transmitters are deployed on corresponding first and second subs61and62that are free to rotate with respect to one another (e.g., in an embodiment in which a drilling motor65is deployed therebetween). As in tool embodiment50, each of the transmitters T1and T2and receivers R1and R1may include a tri-axial antenna arrangement. In the example embodiment depicted the moment of R1zis aligned with the moment of T1z(and the z-axis) while the moments of R1xand R1yare rotationally offset from the moments of T1xand T1yby an offset angle α (e.g., 45 degrees in the depicted embodiment). The moment of R2zis aligned with the moment of T2zwhile the moments of R2xand R2yare rotationally offset from the moments of T2xand T2yby a (e.g., 45 degrees). The disclosed embodiments are, of course, not limited in these regards.

As stated above, the first and second subs61and62may rotate with respect to one another such that the moments of the x- and y-axis transmitting and receiving antennas are misaligned and rotate with respect to one another (i.e., the misalignment angle between the subs varies with time). Using the notation shown onFIG. 2C, at any instant in time, the orientation angle of the x-axis on sub61(the T1xdirection) is θ1with respect to an arbitrary ‘global’ (or wellbore) x-direction. Likewise, at the same instant in time, the orientation angle of the x-axis on sub62(the T2xdirection) is θ2with respect to the global x-direction. It will thus be understood that the moments of the x- and y-transmitting and receiving antennas T1and T2and R1and R2are misaligned by a misalignment angle γ=θ1θ2. It will be understood that θ1and θ2may be referred to as toolface angles of the first and second subs in that they define the rotational orientation of the subs with respect to a global reference direction. Since θ1and θ2are variable with time (owing to the rotation of the subs) and since the subs rotate at different rates the misalignment angle γ also varies with time.

FIG. 3depicts a schematic illustration of a decentered (eccentered) electromagnetic logging tool50deployed in a wellbore40that penetrates an anisotropic formation at a relative dip angle is shown. A wellbore reference frame may be defined by x-, y-, and z-axes (which are fixed relative to the wellbore). A tool reference frame may be defined by x′-, y′-, and z′-axes which are fixed relative to the logging tool. Rotation of the tool in the wellbore causes the x′- and y′-tool axes rotate about the z- and z′-axes with respect to the x- and y-axes of the wellbore. The relative angle θ between the reference frames (e.g., between the x- and x′-axes in the plane orthogonal to the z-axis) is commonly referred to in the art as the toolface angle.

The tool50is shown to be decentered in the wellbore40(having a wellbore diameter d) by a decentering distance decc at a decentering azimuthal angle ψ (in the wellbore reference frame). An apparent decentering azimuth (also referred to as the apparent tool decentering angle AZT) may be defined as the direction of tool decentering in the tool reference frame (e.g., with respect to the x′-axis). The formation is depicted to be anisotropic, having vertical and horizontal conductivities σv and σh at a relative dip angle ϕdipwith respect to the x-axis (i.e., with respect to the wellbore reference frame). An apparent dip azimuth angle is indicated by Φ and represents the relative angle between an orientation marker on the tool (e.g., the x′-axis on the tool which is aligned with the magnetic dipole on the x-axis transmitter) and the direction of the formation's normal vector on the plane orthogonal to the tool's z′-axis. The apparent dip azimuth angle is also referred to herein as the apparent formation azimuth AZF. The conductivity of the drilling fluid is also indicated by σmud.

FIG. 4depicts a flow chart of one disclosed method embodiment100for computing multi-layer formation properties from electromagnetic logging while drilling measurements. An electromagnetic measurement tool (e.g., one of the measurement tools depicted onFIGS. 2B and 2C) is deployed in and rotated in a subterranean wellbore at102(e.g., while drilling the wellbore). Electromagnetic measurements are acquired at104(e.g., via firing the transmitters and receiving the corresponding electromagnetic waves at the receiving antennas) while the tool is rotating. The electromagnetic measurements are processed at106using a first inversion to compute borehole effects (e.g., the effect of the logging tool and the wellbore on the measurements). The borehole effects computed in106are then removed from the electromagnetic measurements acquired at104(e.g., via subtraction) to compute borehole corrected measurements at108. The borehole corrected measurements are then processed with a second inversion at110to compute multi-layer anisotropic formation properties.

With continued reference toFIG. 3, and as known to those of ordinary skill in the art, a time varying electric current (an alternating current) in a transmitting antenna produces a corresponding time varying magnetic field in the local environment (e.g., the tool collar and the formation). The magnetic field in turn induces electrical currents (eddy currents) in the conductive formation. These eddy currents further produce secondary magnetic fields which may produce a voltage response in a receiving antenna (these voltage responses are referred to herein as Vijklwhere i and j represent the transmitter and receiver station and k and l represent the axial orientation of the transmitter and receiver moments).

An electromagnetic logging tool including triaxial transmitters and triaxial receivers (e.g., as depicted onFIGS. 2B and 2C) can obtain a complete set of trans-impedance coupling voltages Vijklusing six antenna firing sub cycles. An example firing sequence for a tool having two transmitter and two receiver stations (i, j=1, 2) is shown below in Table 1.

The transmitting antennas (T1x, T1y, T1z, T2x, T2y, and T2z) may be fired sequentially with each of the six receiving antenna simultaneously measuring a corresponding voltage response (listed in Table 1). In one example embodiment, each transmitter firing is 4 millisecond (0.004 second) followed by a 16 millisecond quiet time for a firing sub-cycle interval of 20 milliseconds. In this example it takes 120 milliseconds (6×20) to acquire a complete set of trans-impedance coupling voltages Vijkl(36 voltage measurements in total as indicated in the table).

One can acquire substantially any suitable number N of sets of the complete 9-component trans-impedance coupling voltages Vijklby repeating the firing sequence N times. It will be understood that using a larger number N of measurement sets tends to give more accurate results with a lower standard deviation due to the random motion of the tool in the borehole, however, a larger number N requires a longer data acquisition time. A maximum value for N can be computed, for example, based on a given drilling speed and a given depth sampling interval. For example only, a hypothetical drilling speed one 100 feet per hour and a six inch depth interval yields 18 seconds (3600/200) for data acquisition per depth interval. Assuming that it takes 120 milliseconds to acquire a single set of trans-impedance voltages yields a maximum value of N=150 (in this example).

Since the logging tool is generally rotating during data acquisition (e.g., at 120 rpm), the toolface angle θ varies with each firing sub-cycle. The toolface angle θ may therefore be measured at each sub-cycle (i.e., with each transmitter antenna firing). The multiple sets of voltage measurements may be referred to herein as Vijklm, m=1, 2, . . . , N and acquired at corresponding toolface angles θijklm. Note that subsets of the coupling voltages have the same toolface angle. For example, V11xxm, V11xym, V11xzm, V12xxm, V12xym, and V12xxmare acquired simultaneously while firing the T1xtransmitter such that θ11xxm, θ11xym, θ11xzm, θ12xxm, θ12xym, and θ12xzmhave the same value (are equal to one another).

In general, the receiving antenna voltages are measured while the tool rotates in the borehole. The measured voltages may be expressed mathematically in terms of their harmonic voltage coefficients, for example, as follows thereby enabling the harmonic voltage coefficients to be obtained:
Vijkl=VDC_ijkl+VFHC_ijklcos(θ)+VFHS_ijklsin(θ)+VSHC_ijklcos(2θ)+VSHS_ijklsin(2θ)  (1)

where VDC_ijklrepresents a DC voltage coefficient, VFHC_ijkland VFHS_ijklrepresent first order harmonic cosine and first order harmonic sine voltage coefficients (also referred to herein as first harmonic cosine and first harmonic sine voltage coefficients), and VSHC_ijkland VSHS_ijklrepresent second order harmonic cosine and second order harmonic sine voltage coefficients (also referred to herein as second harmonic cosine and second harmonic sine voltage coefficients) of the ij transmitter receiver couplings.

The coefficients in Equation 1 can be obtained using a least square curve fitting algorithm from a collection of voltage tool data points VTijklmm=1, 2, . . . N, sampled randomly at θijklm(where VT indicates voltages measured at the tool as opposed to modeled or otherwise computed). This least square fitting process offers an advantage of ensample averaging over the collection of data records used in the fitting to reduce the effect of random tool motion during data acquisition and intrinsic instrumentation random noise.

One example of a least square algorithm to obtain the coefficients is as follows: The array of harmonic voltage coefficients Cijkl(where Cijkl=[VDC_ijkl, VFHC_ijkl, VFHS_ijkl, VSHC_ijkl, VSHS_ijkl]) includes the five coefficients defined in Equation 1. A voltage array Zm=Z(θm, Cijkl), which is a function of the toolface angle and the coefficient array, represents the measured sets of trans-impedance coupling voltages Vijklmdescribed above. Taking the derivative of Z(θm, Cijkl) with respect to Cijklyields the following Gramian matrix:

The least square solution for the coefficient Cijklmay then be given in matrix notation as follows:
Cijkl=(ΓΓT)−1ΓΓ

Co-pending, commonly assigned, and commonly invented U.S. Provisional Patent Application No. 62/222,452, which is incorporated by reference herein in its entirety, discloses a method for computing the mean position of an electromagnetic logging tool in a wellbore undergoing lateral vibrations. The '452 application further demonstrates that the five harmonic voltage coefficients computed above in Cijklclosely match the coefficients that may be obtained with the tool is stationary and located at the mean position. The methods disclosed in the '452 application thus effectively overcome the difficulty of dealing with a tool having a changing and unknown position in the wellbore by solving for the mean tool position and then inverting for formation properties using a model in which the tool is located at the mean position.

Upon obtaining the harmonic voltage coefficients, ratios of the DC xx and yy voltage measurements or the second harmonic xx and yy voltage measurements may optionally be computed and allow a gain ratio of the x to y transmitter and gain ratio of the x to y receiver to be obtained. The voltage measurements may also be rotated mathematically to simulate rotation of the x and y antennas in the R1and R2receivers and the T2transmitter such that they are rotationally aligned with the x and y antennas in the T1transmitter. Such rotation removes the effect of the offset angle α and misalignment angle γ on the measurements. Such computations are disclosed, for example, in U.S. patent application Ser. No. 14/549,396 which is incorporated by reference herein in its entirety.

where theVijklquantities in Equations 6-14 are computed from the harmonic voltage coefficients and may be related to trans-impedance coupling voltages VTijklmeasured by a stationary tool located at the above described mean position in the borehole.

The quantities in Equations 6-14 contain only x and z transmitter and receiver gains. These gains may be canceled out via computing various ones of the following ratios. For example, the following term by term (TBT) compensation operators may be defined for any measurement X obtained between transmitter i and receiver j, for example, as follows:

where Xij, Xji, Xii, and Xjjmay include the measurement terms defined above with respect to Equations 6-14 obtained using the i and j transmitter and receiver (e.g., the transmitters and receivers depicted onFIGS. 2A, 2B, and 2C).

Various gain compensated quantities may be computed following the form of Equation 15. For example, using the TBT operator in Equation 15:

where CVklis a complex compensated quantity and CRVkl, CIVkl, CAVkl, and CPVklare real compensated quantities. The operators real(⋅), imag(⋅), abs(⋅), and p ase(⋅) return the real component, the imaginary component, the absolute value, and the phase angle of the complex quantityVijkl.

As disclosed in more detail in U.S. Provisional Patent Application Ser. No. 62/250,651, which is incorporated by reference herein in its entirety, many of the CVklquantities are azimuthally invariant (e.g., invariant with respect to the formation dip azimuth and/or the decentering azimuth). For example, for a tool having a centered mean position, the off diagonal terms of CVkl, namely CVxz, CVzx, CVyz, CVzy, CVxy, and CVyx, are azimuthally invariant. The following quantities are also azimuthally invariant for a tool having a centered mean position:

where CXXPYY represents a complex measurement quantity related to an xx plus yy coupling component and CXXMYY represents a complex measurement quantity related to an xx minus yy coupling component.

The following measurement quantities are also azimuthally invariant when the mean tool position is off center (decentered):
CRXXPYY=TBT(real(xxpyyij))
CIXXPYY=TBT(imag(xxpyyij))
CRXYMYX=TBT(real(xymyxij))
CIXYMYX=TBT(imag(xymyxij))
CRZ=TBT(real(Vijzz))  (18)

where CRXXPYY and CIXXPYY represent real and imaginary gain compensated quantities related to the xx plus yy coupling, CRXYMYX and CIXYMYX represent real and imaginary gain compensated quantities related to the xy minus yx coupling, and CRZ represents a real gain compensated quantity related to the zz coupling.

Gain compensated quantities that preserve the azimuthal properties may also be computed as disclosed in co-pending, commonly assigned U.S. patent application Ser. No. 14/853,322, which is incorporated by reference herein in its entirety. The following gain compensated quantities preserve the azimuthal properties of the original voltage measurements:

where CZX, CZY, CXZ, and CYZ represent gain compensated quantities that have the characteristics of zx, zy, xz, and yz tensor couplings. These quantities may be further combined to compute gain compensated symmetrized (SX and SY) and anti-symmetrized measurement quantities (AX and AY), as follows:
SX=CZXCXZ
AX=CZX+CXZ
SY=CZYCYZ
AY=CZY+CYZ(20)

It will be understood that the disclosed embodiments are not limited to using the above defined gain compensated measurement quantities. Other suitable gain compensated measurement quantities are disclosed in commonly assigned, co-pending U.S. patent application Ser. Nos. 14/285,581; 14/285,588; 14/339,959; 14/325,797; 14/549,396; and U.S. Provisional Patent Application Ser. No. 62/250,662, each of which is incorporated by reference herein in its entirety.

FIG. 5depicts a flow chart of another disclosed method embodiment200for computing multi-layer formation properties from electromagnetic logging while drilling data. Electromagnetic voltage measurements are acquired at202(e.g., while rotating an electromagnetic logging while drilling tool in a wellbore). The measurements are processed at204using a first inversion in which the inversion model includes tool details and a uniform (single layer) anisotropic formation to output various formation properties, for example, including horizontal and vertical resistivities, a dip angle, a formation dip azimuth, and various wellbore properties, for example, including a mean decentering distance and a mean decentering azimuth of the tool in the wellbore, a resistivity of the drilling fluid, and the hole diameter. A first set of modeled measurements (e.g., modeled voltage measurements) is computed at206using the model (including the tool with finite size antenna coil and conducting mandrel in a borehole through an anisotropic formation) and the inversion output from204. A second set of modeled measurements (e.g., modeled voltage measurements) is computed at208using a point dipole model and the inversion output from204. Differences between the first and second sets of modeled measurements are computed at210to obtain the effects of the borehole (the borehole effects) on the measurements. These borehole effects are subtracted from the acquired measurements at212to obtain borehole corrected measurements. The borehole corrected measurements are then processed using a 1D inversion that includes a point dipole model and a multi-layer anisotropic formation to compute the multi-layer formation properties (e.g., the vertical and horizontal resistivities and a thickness of each of the multiple layers).

FIG. 6depicts a graphical representation of a suitable 1D inversion model for use in the flow chart shown onFIG. 5and the block diagram shown onFIG. 7. The 1D inversion model includes no borehole. The transmitters and receivers are represented mathematically as point dipoles located along a longitudinal tool axis. A multi-layer formation having a fixed depth w includes Nlayer anisotropic formation layers. Each formation layer is defined as having a horizontal and vertical resistivities Rhqand Rvq(q=1, 2, . . . , Nlayer). The bed boundary locations are defined along the length w at locations zq(q=1, 2, . . . , Nlayer) thereby defining the thicknesses of each layer (the thickness of the first layer is z1z0, the thickness of the second layer is z2z1, the thickness of the third layer is z3z2, and so on). The multi-layer formation may be inclined with respect to the tool axis as indicated by formation dip angle DIP and a formation dip azimuth AZF. The formation dip angle and dip azimuth are constant across all layers. The formation resistivities Rhqand Rvq, the bed boundary locations zq, and the constant (or average) formation dip angle DIP and dip azimuth AZF are unknown and may be computed using the 1D inversion.

FIG. 7depicts a block diagram of another disclosed method embodiment300for inverting electromagnetic logging while drilling data to compute multi-layer formation properties. The depicted block diagram may be thought of as including a two-step workflow. In step 1, depicted at305, the electromagnetic logging data (obtained at302) is processed to compute and remove borehole effects. In step 2, depicted at315, the borehole corrected data computed in step 1 is processed using a 1D inversion to compute the multi-layer formation properties.

Electromagnetic voltage measurements are acquired at a computer processor at302(e.g., from an electromagnetic logging while drilling tool that acquired the voltage measurements while rotating in a wellbore as described above). For example, the voltage measurements may be obtained from the logging tool memory upon tripping the tool out of the wellbore or via a high bandwidth communication link (such as wired drill pipe) between the logging tool and a surface processor. A large number of voltage measurements are generally acquired. These measurements are represented in302by VTijklmp, where VT represent tool voltage measurements, i and j represent the transmitter and receiver stations (e.g., 1, 2), k and l represent the transmitting antenna and receiving antenna moment axes (e.g., x, y, z), m represents the voltage set (sub-cycle) at a single depth such that m=1, 2, . . . , N, and p represents the number of depth cycles such that p=1, 2, . . . , ndept. Measured or estimated drilling fluid resistivity values RMpand borehole diameter values HDp(where p is as defined above) are also acquired at304.

At306(in step 1), the acquired tool voltage measurements are processed using an inversion model that includes tool and borehole details as well as a uniform anisotropic formation to compute various tool, borehole, and formation properties. For example, the inversion may be used to compute the horizontal and vertical resistivities of the uniform anisotropic formation Rhand Rv, the mean position of the logging tool in the borehole deccmeanand Ψmean, the formation dip and dip azimuth DIP and AZF, the resistivity of the drilling fluid RM, and the borehole diameter HD.

In block308, modeled trans-impedance voltage responses VMijklmare computed using the forward model and the outputs from block306, where VM represents the modelled voltages. The modeled voltage response are further processed to compute a plurality Nc of compensated measurement quantities CMn(where CM represents the modeled compensated measurement quantities and n=1, 2, . . . , Nc). These quantities may be computed, for example, as described above in Equations 15-20 and/or in the various cited references.

In block310modeled trans-impedance voltage responses VMNBHijklmare computed using the formation property outputs (Rh, Rv, DIP, and AZF) from block306and a point dipole model with no borehole and point dipole transmitting and receiving antennas. VMNBH represents the modelled voltages in the absence of a borehole. The modeled voltage responses are further processed to compute a plurality Nc of compensated measurement quantities CMNBHn(where CMNBH represents the modeled compensated measurement quantities in the absence of a borehole and n=1, 2, . . . , Nc). These quantities may be computed, for example, as described above in Equations 15-20 and/or in the various cited references.

At block312the borehole effects on the modeled compensated measurement quantities may be computed, for example, as follows:
ΔCMn=CMnCMNBHn, n=1,2, . . . ,Nc(21)

where ΔCMnrepresents the borehole effects on each of the n=1, 2, . . . , Nc compensated measurement quantities.

At block314, compensated measurements CTnare computed from the tool voltage measurements VTijklmthat were received at306. These compensated quantities are further processed in combination with the borehole effects ΔCMncomputed in312to compute borehole corrected compensated quantities CBHCn, for example, as follows:
CBHCn=CTnΔCMn, n=1,2, . . . ,Nc(22)

With respect to blocks306,308,310,312, and314(which collectively make up step 1 as depicted at305), it will be understood that these blocks are repeated at each depth interval (at each depth p) such that the above described borehole and formation parameters, modeled voltages, and compensated quantities are computed independently at each depth (i.e., at p=1, p=2, and so on). Thus, each of these quantities may also be represented as including a subscript p, for example, CBHCnprepresenting the n=1, 2, . . . , Nc borehole corrected compensated measurement quantities at each of the p=1, 2, . . . , ndept depth intervals.

As stated above, step 2 (depicted at315) includes processing the borehole corrected measurements using a 1D inversion to compute the multi-layer formation properties. In particular, the compensated borehole corrected measurements CBHCnpare processed at316using a 1D inversion in which the forward model includes point dipole transmitters and receivers without a tool body and without a borehole in a multi-layered anisotropic formation. The outputs318include the horizontal and vertical resistivities of each layer Rhqand Rvq, the bed boundary locations zq(q=1, 2, . . . , Nlayer), and the constant (or average) formation dip angle DIPaveand dip azimuth AZFave.

One advantageous feature of method300is that the borehole correction process in step 1 is configured to match the inversion engine in step 2 such that the correction acts as a bridge between the inversion model in block306and the layered formation model in block316. In particular, matching point dipole models are used in blocks310and316. In a uniform anisotropic formation the compensated tool measurements CTnmay be about equal to the modeled compensated measurements CMnsuch that CTn˜CMn. Thus, from Equations 21 and 22 the compensated borehole corrected measurements CBHCnmay be about equal to the point dipole modeled compensated measurements CMNBHnsuch that CBHCn=CTnΔCMn˜CMNBHn.

Another advantageous feature of method300is that the inversion engines employ fully gain compensated measurement quantities (e.g., as opposed to trans-impedance voltages or an apparent conductivity matrix). The borehole correction process described above with respect toFIG. 7involves three distinct responses; (i) the actual tool measurements, (ii) model responses from a finite size antenna coil in a borehole, and (iii) model responses from a point dipole tool with no borehole in a layered formation. The magnitudes of these depend on the gains of the antennas in the tool and the models. In practice, the gains depend on the tool operating pressure and temperature (as well as other factors) and thus can change during a drilling operation. The compensated measurements disclosed herein are advantageously independent of gain (and gain variations) and therefore may provide accurate borehole effect measurements and corrections as well as accurate formation properties upon inversion.

FIG. 8depicts a block diagram that provides further detail regarding element306of one embodiment ofFIG. 7. Block322may be referred to as the forward model engine of the first inversion. It computes the model responses (model compensated quantities) for a given set of borehole/formation parameters, namely, in this example: Rh, Rv, deccmean, Ψmean, DIP, AZF, RM, and HD. Of this list of parameters, some may be treated as known and thus not inverted (such as Ψmean, AZF, RM, and HD). The remaining parameters (e.g., Rh, Rv, decc, and DIP) are referred to as free parameters to be inverted. As described above, the model includes tool and borehole details as well as a uniform anisotropic formation. Finite element code may be used to generate the model responses (although such code tends to be too slow for real time while drilling applications). Alternatively, finite element code may be used to pre-compute a multi-dimensional response table in which each dimension corresponds to a borehole/formation parameter. The response table may then be stored in a computer so that the model responses at arbitrary borehole/formation parameter values may be quickly computed through interpolation. The disclosed embodiments are not limited in these regards.

Compensated measurement quantities are computed in block324using the acquired trans-impedance voltage measurements (at302inFIG. 7). It will be understood that to process a proper inversion the number of linearly independent compensated measurements generated in block324generally needs to exceed the number of free parameters in block322. Substantially any number of the disclosed compensated measurement quantities may be computed in block324.

The model quantities computed in322may be compared with the tool quantities computed in324using a cost function in block326. The cost function may be configured to have a decreasing value as the model measurements approach the tool measurements. The cost function may be further configured to equal zero when the model measurements exactly match the corresponding tool measurements. Substantially any suitable cost function may be utilized, for example, including a quadratic form of the difference between the tool measurements and model responses. When the cost function is less than a predetermined threshold at328the present values of the borehole and formation parameters Rh, Rv, deccmean, Ψmean, DIP, AZF, RM, and HD are output at334. Otherwise, the free parameters (e.g., Rh, Rv, deccmean, and DIP) may be adjusted at332and322the inversion loop repeated beginning at322as depicted.

FIG. 9depicts a block diagram that provides further detail regarding element316of one embodiment ofFIG. 7. Blocks342and344describe inputs to the 1D inversion. Block342includes the borehole corrected compensated measurements CBHCnpfrom block314ofFIG. 7. Note that there are n=1, 2, . . . , Nc borehole corrected compensated measurement quantities at each of the p=1, 2, . . . , ndept depths. Block344receives Rhp, Rvp, DIPpand AZFpvalues from block306ofFIG. 7(note that these measurements are received for each depth interval and thus include the subscript p).

Block346processes Rhp, Rvp, DIPpand AZFpvalues from block344to compute initial estimated parameters for the 1D inversion. With reference again toFIG. 6, the estimated parameters may include the horizontal and vertical resistivities of each of the layers, Rhqand Rvq, the locations (depths) of the bed boundaries zq, and the average dip angle and dip azimuth DIPaveand AZFave.

In one example embodiment a first method for establishing the initial estimates is employed when the average dip angle DIPaveis less than a threshold (a low angle well condition) and a second method for establishing the initial estimates is employed when the DIPaveis greater than the threshold (a high angle well condition). The average dip angle DIPavemay be computed, for example as follows:

Substantially any threshold may be used for distinguishing between low and high angle well conditions. For example only, the threshold may be 75 degrees such that low angle conditions apply when DIPave≤75° and high angle conditions apply when DIPave>75°. Generally, during low angle conditions the logging tool traverses multiple formation layers and therefore multiple boundaries (e.g., as depicted onFIG. 6). During high angle conditions the logging tool generally traverses only a few (or no) layers/boundaries. As DIPaveapproaches 90 degrees, the logging tool may reside in a single layer.

FIG. 10Adepicts an example well path402versus depth w through multiple formation layers and boundaries (where w represents depth on a total vertical depth TVD scale). This example is representative of a low angle operation in which the well path traverses several layers. A log404of the horizontal resistivity Rhversus w is also depicted.

In this low angle example, the bed boundary depths zqmay be estimated by solving for the locations (depths) of the inflection points in the Rhlog404. This may be accomplished by solving for the depths w at which the second derivative of the Rhlog404with respect to depth w is equal to zero, for example, as follows:

Note that the number of layers selected in the 1D inversion is indeterminate (not fixed), but may be selected based on the measured data, for example, based on the number of zeros in the second derivative of the Rhlog. The initial estimates for Rhqand Rvqmay be obtained by computing the inverse of the averaged horizontal and vertical conductivities in each of the layers, for example, as follows:

where s=1, 2, . . . , Nqrepresent the depth sample points with the q-th layer. Rhqis depicted onFIG. 10Aas the vertical dotted lines (at406). The average formation dip azimuth AZFavemay be computed, for example, as follows:

where a tan 2(⋅) represents the four quadrant inverse tangent function and MOD(⋅) represents the modulo operation which gives the remainder after division.

FIG. 10Bdepicts an example well path412versus depth w (where w represents depth on a total vertical depth TVD scale). This example is representative of a high angle operation in which the well path remains in a single layer. A log414of the horizontal resistivity Rhversus w is also depicted. In this example, the variation in Rhis primarily due to the variation in the distance of the well path412from the nearest bed boundary.

In a high angle operation it may be advantageous to utilize a 1D inversion model employing a fixed number of layers (since the wellbore likely crosses few, if any, boundaries). For example, the 1D inversion model may arbitrarily be selected to include 2, 3, or 4 layers. The initial estimates of the horizontal and vertical resistivities of each layer may be obtained, for example, by simply averaging the values over the entire depth cycle as follows:

Moreover, the bed boundary locations (depths) may be selected arbitrarily, for example, by dividing the depth w into equal intervals. The average resistivity Rh,avemay also be used to determine a maximum range of zq. A higher Rh,avevalue implies that the tool is in a high resistivity environment with a larger range of investigation, i.e., a larger range for zq. Conversely, a lower Rh,aveimplies that a smaller range of zqmay be appropriate.

With reference again toFIG. 9, block348, which depicts the forward 1D model engine, is configured to compute model compensated measurement quantities CMnpbased upon the current values of Rhq, Rvq, zq, DIPave, and AZFave. Various analytical methods are known for making such computations, for example, as disclosed by Anderson in “Modeling and Inversion Methods for the Interpretation of Resistivity Logging Tool Response”, p. 285-293, DUP Science, ISBN 90-407-2231-5.

Block350compares the modeled compensated measurements CMnpwith the borehole corrected compensated measurements CBHCnpinput in342using a cost function. The cost function may be configured to have a decreasing value when the model measurements more closely match the corrected tool measurements. The cost function may be further configured to equal zero when the model measurements exactly match the corresponding tool measurements. Example cost functions may be represented mathematically, for example, as follows:

where C represents the cost function, Nc represents the number of compensated measurements, ndept represent the number of depth intervals and wnprepresent weighing coefficients of the n-th compensated measurement at the p-th depth interval. The cost function may further optionally include regularization terms such that it generates smoother outputs for noisy input data condition, although the disclosure is by no means limited in this regard.

When the cost function is less than a predetermined threshold at352the present values of the multi-layer formation parameters Rhq, Rvq, zq, DIPave, and AZFaveare output at356. Various inversion parameters may also be output at356. Otherwise, these parameters are adjusted at354and then used to compute new modeled compensated measurement quantities at348.

It will be understood that blocks348,350,352, and354form an inversion loop. The number of times the algorithm loops through block352, Nloop, may be recorded and used for determining loop termination in the event that there is no convergence. Loop exiting criteria may include, for example, (i) the number of iterations Nloopexceeds a threshold, (ii) the cost function is less than a threshold indicating convergence, and (iii) a difference between successive cost function computations is less than a threshold indicating convergence.

It will be understood that the various methods disclosed herein for inverting electromagnetic logging data may be implemented on a on a computer, for example, including one or more microprocessors and electronic memory. The computer may be configured to acquire the logging data from an electromagnetic logging tool, for example, after the tool has been removed from the subterranean wellbore (e.g., via a communications link). Alternatively, the computer may be configured to acquire logging data from the electromagnetic logging tool in substantially real time while the tool is making measurements (e.g., via high bandwidth wired drill pipe telemetry system). The computer may be further configured to execute the various method embodiments discloses herein, for example, inFIGS. 4, 5, 7, 8, and 9, and process one or more of the disclosed mathematical equations.

Although methods for inverting electromagnetic logging measurements have been described in detail, it should be understood that various changes, substitutions and alternations can be made herein without departing from the spirit and scope of the disclosure as defined by the appended claims.