Abstract:
A system and method for estimating porosity distribution in a region of interest of a geologic formation from a resistivity image log representative of the geologic formation is disclosed. A normalization factor representative of a rock matrix based on a first resistivity value and an image point factor based on a second resistivity value are calculated and compared to identify points in the resistivity image log that correspond to the secondary porosity. The normalization factor and image point factor are recalculated based on a different first resistivity value and a different second resistivity value as necessary to identify additional points in the resistivity image log that correspond to the secondary porosity until a termination criterion is met. The method may further include a porosity calibration operation and one or more artifact corrections.

Description:
RELATED APPLICATIONS 
       [0001]    This application is a divisional application of U.S. patent application Ser. No. 13/945,690, filed Jul. 18, 2013 entitled “System and Method For Estimating Porosity Distribution in Subterranean Reservoirs” which is hereby incorporated by reference in its entirety. 
     
    
     STATEMENT REGARDING FEDERALLY SPONSORED RESEARCH OR DEVELOPMENT 
       [0002]    Not applicable. 
       FIELD OF THE INVENTION 
       [0003]    The present invention relates generally to methods and systems for processing well logs and, in particular, methods and systems for estimating porosity distribution, including secondary porosity, in subterranean reservoirs. 
       BACKGROUND OF THE INVENTION 
       [0004]    Hydrocarbon exploration and production is aided by understanding the porosity distribution in subterranean reservoirs. In some geologic formations, particularly carbonate formations, the porosity distribution may be significantly non-uniform on one or more scales from core plug scale to interwell distances. 
         [0005]    Porosity may be estimated by inspection of core samples or evaluation of well logs. However, these approaches have difficulties when pore space structure changes on a scale shorter than the spacing of the core measurements or shorter than the log sensitivity. In the presence of very large pores (vugs, for example), core samples may not be large enough to capture a single large vug occurrence nor a representative distribution of vugs needed to characterize fluid flow on a well-scale. These problems are common in many carbonate fields. 
       SUMMARY OF THE INVENTION 
       [0006]    Described herein are implementations of various approaches for a computer-implemented method for estimating porosity distribution in regions of interest of geologic formations. 
         [0007]    A computer-implemented method for estimating porosity distribution in a region of interest of a geologic formation from a resistivity image log representative of the geologic formation including calculating a normalization factor representative of a rock matrix based on a first resistivity value; calculating a image point factor based on a second resistivity value; comparing the image point factor and the normalization factor to identify points in the resistivity image log that correspond to the secondary porosity; recalculating the normalization factor and the image point factor based on a different first resistivity value and a different second resistivity value; re-comparing the recalculated normalization factor and the recalculated image point factor to identify additional points in the resistivity image log that correspond to the secondary porosity; and repeating the recalculating and re-comparing steps until a termination criterion is met is disclosed. The method may further include a porosity calibration operation and one or more artifact corrections. 
         [0008]    In another embodiment, a computer system including a data source or storage device, at least one computer processor and a user interface to implement the method for estimating porosity distribution in a region of interest of a geologic formation is disclosed. 
         [0009]    In yet another embodiment, an article of manufacture including a computer readable medium having computer readable code on it, the computer readable code being configured to implement a method for estimating porosity distribution in a region of interest of a geologic formation is disclosed. 
         [0010]    The above summary section is provided to introduce a selection of concepts in a simplified form that are further described below in the detailed description section. The summary is not intended to identify key features or essential features of the claimed subject matter, nor is it intended to be used to limit the scope of the claimed subject matter. Furthermore, the claimed subject matter is not limited to implementations that solve any or all disadvantages noted in any part of this disclosure. 
     
    
     
       BRIEF DESCRIPTION OF THE DRAWINGS 
         [0011]    These and other features of the present invention will become better understood with regard to the following description, claims and accompanying drawings where: 
           [0012]      FIG. 1  is a diagram of porosity in a geologic formation; 
           [0013]      FIG. 2  is a flowchart of an embodiment of the invention; 
           [0014]      FIG. 3A  is a diagram of a resistivity imaging tool; 
           [0015]      FIG. 3B  is a diagram of a part of the resistivity imaging tool; 
           [0016]      FIG. 4  illustrates an intermediate step of an embodiment of the invention; 
           [0017]      FIG. 5  illustrates results for porosity distribution and secondary porosity estimate from an embodiment of the invention; and 
           [0018]      FIG. 6  schematically illustrates a system for performing a method in accordance with an embodiment of the invention. 
       
    
    
     DETAILED DESCRIPTION OF THE INVENTION 
       [0019]    The present invention may be described and implemented in the general context of a system and computer methods to be executed by a computer. Such computer-executable instructions may include programs, routines, objects, components, data structures, and computer software technologies that can be used to perform particular tasks and process abstract data types. Software implementations of the present invention may be coded in different languages for application in a variety of computing platforms and environments. It will be appreciated that the scope and underlying principles of the present invention are not limited to any particular computer software technology. 
         [0020]    Moreover, those skilled in the art will appreciate that the present invention may be practiced using any one or combination of hardware and software configurations, including but not limited to a system having single and/or multiple processor computers, hand-held devices, tablet devices, programmable consumer electronics, mini-computers, mainframe computers, and the like. The invention may also be practiced in distributed computing environments where tasks are performed by servers or other processing devices that are linked through one or more data communications networks. In a distributed computing environment, program modules may be located in both local and remote computer storage media including memory storage devices. The present invention may also be practiced as part of a down-hole sensor or measuring device or as part of a laboratory measuring device. 
         [0021]    Also, an article of manufacture for use with a computer processor, such as a CD, pre-recorded disk or other equivalent devices, may include a tangible computer program storage medium and program means recorded thereon for directing the computer processor to facilitate the implementation and practice of the present invention. Such devices and articles of manufacture also fall within the spirit and scope of the present invention. 
         [0022]    Referring now to the drawings, embodiments of the present invention will be described. The invention can be implemented in numerous ways, including, for example, as a system (including a computer processing system), a method (including a computer implemented method), an apparatus, a computer readable medium, a computer program product, a graphical user interface, a web portal, or a data structure tangibly fixed in a computer readable memory. Several embodiments of the present invention are discussed below. The appended drawings illustrate only typical embodiments of the present invention and therefore are not to be considered limiting of its scope and breadth. 
         [0023]    The present invention relates to estimating porosity distribution in a geologic formation, particularly in a carbonate formation with secondary porosity (such as without limitation, vugs, molds or dissolution enhanced fractures). Significantly non-uniform porosity distribution is common in carbonate reservoirs on all lengthscales of routine measurements in oil exploration and production, from core plug scale to interwell distances. An accurate representation of porosity is desirable in building reservoir models for estimation of oil-in-place and recoverable reserves. 
         [0024]    In rocks with simple, uniform pore systems filled with saline water, Archie&#39;s law relates porosity φ to the formation resistivity factor F, as 
         [0000]    
       
         
           
             F 
             = 
             
               1 
               
                 φ 
                 
                   
                       
                   
                    
                   m 
                 
               
             
           
         
       
     
         [0000]    where F is the ratio of resistivity R 0  of fully water-saturated rock and resistivity of formation water R W   
         [0000]    
       
         
           
             F 
             = 
             
               
                 R 
                 0 
               
               
                 R 
                 W 
               
             
           
         
       
     
         [0000]    and m is the lithology exponent or cementation factor. Those skilled in the art are aware that, when hydrocarbons are present, the resistivity index I, which is a ratio of resistivity of the rock containing hydrocarbons R t  to the resistivity of the fully water-saturated rock R 0 , i.e. I=R t /R 0 , is related to the water saturation S w  of the rock as 
         [0000]    
       
         
           
             I 
             = 
             
               1 
               
                 S 
                 w 
                 n 
               
             
           
         
       
     
         [0000]    where n is the saturation exponent. In sandstones, it has been shown that m and n are very nearly equal to 2, and using the two relations the well-known Archie equation is obtained: 
         [0000]    
       
         
           
             
               S 
               w 
               n 
             
             = 
             
               
                 
                   R 
                   W 
                 
                 
                   R 
                   t 
                 
               
               / 
               
                 φ 
                 m 
               
             
           
         
       
     
         [0025]    At least two conditions exist in which the saturation exponent n can be significantly different than 2: in oil-wet reservoirs where oil coats grains and starts blocking pore throats and suppressing electrical conduction, even at low oil saturation, and in rocks with non-uniform pore space. The cementation factor m is related to the tortuosity of current paths, taking a value of 1 in an idealized reservoir where fractures offer a straight conductive path without interaction with granular (inter-grain) porosity. On the other hand, isolated pores have porosity but do not contribute to rock conductivity, effectively raising F and m for the host rock. 
         [0026]    In carbonates, approaches which determine coefficients m and n from core or log measurements have been used. These approaches may have difficulties when the pore-space structure changes on a lengthscale shorter than the spacing of core measurements or shorter than the log sensitivity. In the presence of very large pores (vugs, for example), core samples may not be large enough to capture a single vug occurrence nor a representative distribution of vugs needed to characterize fluid flow on a wellscale. Both occurrences, of variations in pore-space structure occurring on small lengthscales and occurrence of vugs, are very common in many carbonate fields.  FIG. 1  shows a representative diagram of porosity in a carbonate formation  10 . The rock matrix has pores  11 . There may be secondary porosity such as vugs  12 . In some locations, the porosity may be very low, essentially zero, such as the tight rock  13 . When these tight rock regions occupy layers of thickness larger than the resolution of porosity logs, their presence can be found from the total porosity log being very low (for example below 1%). This example is not meant to be limiting as the determination of a region of tight rock may be complicated by artifacts in the porosity log. In the cases of some minerals, tight regions can be identified from core and a minimum resistivity image value in such regions can be inferred from the image and used as a cut-off to identify such tight regions from image in other regions of the well. 
         [0027]    As indicated in method  15  of  FIG. 2 , the present invention obtains resistivity image log(s) (operation  20 ) from which the porosity distribution is estimated. Resistivity image logs have a resolution significantly higher than conventional open-hole resistivity logs and display results of an array of electrode measurements around the borehole as a depth and azimuth-dependent image. An example of a resistivity image tool may be seen in  FIGS. 3A and 3B . 
         [0028]    The probe may be a multi-trace or multi-pad measurement probe. For example,  FIG. 3A  illustrates a probe  100  for use in borehole characterization that includes a generally elongated shaft  120  having at one end a number of outwardly extending members  140 . The outwardly extending members  140  may each include a pad  160  (shown in more detail in  FIG. 3B ) for interrogating a region of a borehole. The illustrated pad  160  includes a plurality of pairs of sensors  200  for monitoring a current which flows as a result of applying an alternating exciting voltage between an electrode located elsewhere on the tool. 
         [0029]    In use, the probe  100  is generally lowered into the borehole to be characterized. Upon reaching an appropriate depth, which may be the bottom of the hole, or a selected intermediate depth, the probe is retrieved and measurements are taken as the probe rises through the material. In many cases, the probe  100  will have four pads  160  so that the hole may be characterized in four regions with distinct azimuths. In another example, the probe  100  may have six pads  160  that characterize regions around six distinct azimuths. The pad  160  may be accompanied by a flap  260  which also has a plurality of pairs of sensors  200 . 
         [0030]    The sensors  200  on the pad  160  and/or the flap  260  measure the electrical current that passes through the geologic formation and is proportional to the formation conductivity. Each sensor  200  will measure the current in its immediate vicinity, thereby measuring the conductivity of the geologic formation directly facing it. Resistivity is inversely proportional to the measured current and the coefficient of proportionality is the same for all electrodes within a homogeneous region. The coefficient of proportionality is well approximated as the same for all electrodes in other cases provided the number of electrodes is large. Measurement results are processed to obtain an array referred to in the art as the raw image. For better contrast in viewing, the image is often processed further. However, to retain relation of the image values to resistivity, only a calibration of the raw image to a conventional shallow-resistivity measurement is needed, without contrast enhancement procedures. 
         [0031]      FIG. 4  shows an example of a raw image log in column  52 , from data recorded by a Schlumberger Fullbore MicroImager (FMI) resistivity imaging tool. The four measurements from the four pad-flaps appear as individual columns at locations in image corresponding to the azimuth of their measurement location in the borehole. Column  50  represents the depth in the borehole (with reference altered for confidentiality). Column  54  shows the calibrated image. Column  56  shows the calibration value as the dashed dark line and the resultant average calibrated image value as the light solid line. The calibration value is also shown in gray-scale in column  57  and the average calibrated image value is shown in column  58 . 
         [0032]    Methods have been developed to model porosity for each point in the image using Archie&#39;s equation and may include using constraints from other porosity information and a conventional resistivity measurement. One method estimates the fraction of porosity which corresponds to regions of secondary porosity by examining the distribution of image-derived porosity and statistically defining cutoff values on porosity. However, this type of analysis is not appropriate for conductive minerals and shale. This method may overestimate secondary porosity in regions where the rock matrix has significant porosity variations with abrupt onset azimuthally (around the borehole). 
         [0033]    Referring again to  FIG. 2 , obtaining the resistivity image log (operation  20 ) may include running the probe in the borehole or receiving the data recorded by the tool, processing it to obtain a raw image and calibrating the raw image. Once the data has been obtained, the present invention assumes that the rock matrix obeys Archie&#39;s law or a similar relation 
         [0000]    
       
         
           
             
               
                 
                   
                     
                       φ 
                       i 
                     
                     = 
                     
                       
                         
                           ( 
                           
                             
                               1 
                               
                                 S 
                                 w 
                                 n 
                               
                             
                              
                             
                               
                                 R 
                                 w 
                               
                               
                                 r 
                                 i 
                               
                             
                           
                           ) 
                         
                         
                           1 
                           / 
                           m 
                         
                       
                       ≡ 
                       
                         C 
                         
                           r 
                           i 
                           
                             1 
                             / 
                             m 
                           
                         
                       
                     
                   
                   , 
                 
               
               
                 
                   ( 
                   1 
                   ) 
                 
               
             
           
         
       
     
         [0000]    where φ i  and r i  are respectively the rock matrix porosity and resistivity in the i-th cell contributing to the i-th electrode of the fine-resolution resistivity measurement, R w  is the resistivity of the saline water in the matrix pores, and C is a constant for constant water saturation S w , salinity and temperature around the borehole at the given reference depth. The rock matrix excludes regions of tight rock where porosity can be assumed to be zero ( FIG. 1 , tight rock  13 ). 
         [0034]    Expressing the total volume of voids in a rock via volume of pores in the matrix and of voids in the region affected by secondary porosity (vugs  12  in  FIG. 1 ), the constant C, which is set for a constant water saturation, salinity and temperature around the borehole at a given reference depth, can be related to a reference porosity log φ: 
         [0000]    
       
         
           
             
               
                 
                   
                     φ 
                     × 
                     
                       V 
                       total 
                     
                   
                   = 
                   
                     
                       
                         V 
                         cell 
                       
                       × 
                       
                         
                           ∑ 
                           
                             i 
                             = 
                             1 
                           
                           
                             N 
                             ma 
                           
                         
                          
                         
                           C 
                           
                             r 
                             i 
                             
                               1 
                               / 
                               m 
                             
                           
                         
                       
                     
                     + 
                     
                       V 
                       sec 
                     
                   
                 
               
               
                 
                   ( 
                   2 
                   ) 
                 
               
             
           
         
       
     
         [0000]    where φ may be determined from core measurements or well logs such as a neutron-density crossplot, V sec  is the volume sensed by the fine-resolution measurement as occupied by secondary porosity, V total  is the volume of the borehole region approximately shaped as a cylindrical shell, which is probed by the tool which provided reference porosity at the given depth and comprises of non-overlapping cells each sensed by one electrode of the fine-resolution resistivity measurement, V cell  is the volume of the region (cell) most directly probed by the resistivity measurement of high resolution (i.e. for a resistivity image tool this is exclusive volume assigned to be responsible for a given pixel in image, though the sensitive volume for the tool is larger). In an embodiment, we assume that the ratio V sec /V total  is well approximated by the ratio of the number of resistivity image pixels occupied by secondary porosity to the total number of resistivity image pixels for the region. This assumption is not to be taken as limiting the scope of the invention, as the ratios under consideration can be assumed to be related with another coefficient of proportionality characteristic of the region or formation. Further, N ma  is the number of cells occupied by matrix, and V total  and V cell  are related via the total number of cells N as 
         [0000]        V   total   =N×V   cell   (3)
 
         [0000]    while V total  and V sec  are related as 
         [0000]    
       
         
           
             
               
                 
                   
                     
                       N 
                       ma 
                     
                     N 
                   
                   = 
                   
                     
                       1 
                       - 
                       
                         
                           N 
                           sec 
                         
                         N 
                       
                       - 
                       
                         
                           V 
                           tight 
                         
                         
                           V 
                           total 
                         
                       
                     
                     ≡ 
                     
                       1 
                       - 
                       v 
                       - 
                       
                         
                           V 
                           tight 
                         
                         
                           V 
                           total 
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   4 
                   ) 
                 
               
             
           
         
       
     
         [0000]    where V tight  is the volume of the region occupied by tight rock of essentially no porosity and N sec  is the number of cells that are sensed as occupied by secondary porosity by a fine-resolution resistivity measurement. We call ν=N sec /N the fractional volume sensed as secondary porosity. In a measurement of a resolution much finer than the characteristic size of macropores involved in secondary porosity, N sec =V sec /V cell  and N sec /N=V total . However, for a measurement of resolution of the same order as the size of macropores sensed, the measurement may detect the presence of a void in a region where porosity is significantly less than 100% for a significant portion the total number of pixels which are sensing that void. 
         [0035]    Dividing the left and right hand side of the equation (2) with V total , and using equations (3) and (4) one obtains 
         [0000]    
       
         
           
             
               
                 
                   φ 
                   = 
                   
                     
                       
                         ( 
                         
                           1 
                           - 
                           v 
                           - 
                           
                             
                               V 
                               tight 
                             
                             
                               V 
                               total 
                             
                           
                         
                         ) 
                       
                       × 
                       
                         1 
                         
                           N 
                           ma 
                         
                       
                       × 
                       
                         
                           ∑ 
                           
                             i 
                             = 
                             1 
                           
                           
                             N 
                             ma 
                           
                         
                          
                         
                           C 
                           
                             r 
                             i 
                             
                               1 
                                
                               
                                 / 
                               
                                
                               m 
                             
                           
                         
                       
                     
                     + 
                     
                       φ 
                       v 
                     
                   
                 
               
               
                 
                   ( 
                   5 
                   ) 
                 
               
             
           
         
       
     
         [0000]    Here we call φ ν  secondary porosity as it is the ratio of volume occupied by secondary porosity to the total volume. This can be related to the fractional volume ν as φ ν =ν×φ high  where φ high  is the average porosity assigned to the region in which the fine-resolution resistivity measurements can sense the void. In a measurement of resolution much finer than the size of regions occupied by secondary porosity, φ high =1 but for resolution of the same order as for example one of dimensions of vugs sensed, this parameter is about 50%. These numbers are examples of common parameter values but are not meant to be limiting; other average porosity values are possible and fall in the scope of the present invention. 
         [0036]    From equations (1) and (5), the porosity of a rock matrix cell can be expressed as 
         [0000]    
       
         
           
             
               
                 
                   
                     
                       φ 
                       i 
                     
                     = 
                     
                       
                         
                           φ 
                           - 
                           
                             φ 
                             v 
                           
                         
                         
                           1 
                           - 
                           v 
                           - 
                           
                             
                               V 
                               tight 
                             
                             
                               V 
                               total 
                             
                           
                         
                       
                       × 
                       
                         
                           1 
                           
                             r 
                             i 
                             
                               1 
                                
                               
                                 / 
                               
                                
                               m 
                             
                           
                         
                         
                           
                             〈 
                             
                               1 
                               
                                 r 
                                 
                                   1 
                                    
                                   
                                     / 
                                   
                                    
                                   m 
                                 
                               
                             
                             〉 
                           
                           ma 
                         
                       
                     
                   
                   , 
                 
               
               
                 
                   ( 
                   6 
                   ) 
                 
               
             
           
         
       
     
         [0000]    where          r 1/m             ma  is the average of inverse mth root of the fine-resolution resistivity over the matrix cells: 
         [0000]    
       
         
           
             
               
                 
                   
                     
                       〈 
                       
                         1 
                         
                           r 
                           
                             1 
                              
                             
                               / 
                             
                              
                             m 
                           
                         
                       
                       〉 
                     
                     ma 
                   
                   = 
                   
                     
                       1 
                       
                         N 
                         ma 
                       
                     
                      
                     
                       
                         ∑ 
                         
                           j 
                           = 
                           1 
                         
                         
                           N 
                           ma 
                         
                       
                        
                       
                         1 
                         
                           r 
                           j 
                           
                             1 
                              
                             
                               / 
                             
                              
                             m 
                           
                         
                       
                     
                   
                 
               
               
                 
                   ( 
                   7 
                   ) 
                 
               
             
           
         
       
     
         [0000]    The equation (6) is not applicable in cases where 
         [0000]    
       
         
           
             
               
                 
                   
                     
                       
                         
                           ( 
                           
                             φ 
                             - 
                             
                               φ 
                               v 
                             
                           
                           ) 
                         
                         × 
                         
                           1 
                           
                             r 
                             i 
                             
                               1 
                                
                               
                                 / 
                               
                                
                               m 
                             
                           
                         
                       
                       
                         
                           ( 
                           
                             1 
                             - 
                             v 
                             - 
                             
                               
                                 V 
                                 tight 
                               
                               
                                 V 
                                 total 
                               
                             
                           
                           ) 
                         
                         × 
                         
                           
                             〈 
                             
                               1 
                               
                                 r 
                                 
                                   1 
                                    
                                   
                                     / 
                                   
                                    
                                   m 
                                 
                               
                             
                             〉 
                           
                           ma 
                         
                       
                     
                     &gt; 
                     1 
                   
                   , 
                 
               
               
                 
                   ( 
                   8 
                   ) 
                 
               
             
           
         
       
     
         [0000]    which can occur with a high proportion occupied by tight regions or for pixels in highly conductive regions (where the value of r i  is too low relative to its peers considered to be in the matrix). It is common to impose cut-offs on resistivity or conductivity to delineate regions of tight rock or regions of secondary porosity. The practical aspect of the restriction in the inequality (8) is that it can be used iteratively to find the maximum resistivity cut-off for secondary porosity voids when the region of tight rock is not present or has already been delineated, e.g. with a cut-off r tight  imposed as a minimum resistivity for a point to belong to the tight rock region. Starting from the lowest resistivity r i  in a region for which porosity distribution is calculated, the numerator in the inequality (8), which can be called the image point or cell factor, is found. Referring again to  FIG. 2 , the image point factor is calculated at operation  24  as 
         [0000]    
       
         
           
             
               
                 ( 
                 
                   φ 
                   - 
                   
                     φ 
                     v 
                   
                 
                 ) 
               
               × 
               
                 1 
                 
                   r 
                   i 
                   
                     1 
                      
                     
                       / 
                     
                      
                     m 
                   
                 
               
             
             , 
           
         
       
     
         [0000]    with the secondary porosity φ ν =0 corresponding to no cells yet assigned to secondary porosity in the first iteration. Similarly the denominator in the inequality (8), which can be called the normalizing factor, can be calculated, as shown in  FIG. 2  operation  22 , with all cells (except for tight regions) belonging to the matrix, i.e. minimum resistivity taken into the average set to the lowest resistivity image value r min  in the region: 
         [0000]    
       
         
           
             
               ( 
               
                 1 
                 - 
                 v 
                 - 
                 
                   
                     V 
                     tight 
                   
                   
                     V 
                     total 
                   
                 
               
               ) 
             
             × 
             
               
                 
                   〈 
                   
                     1 
                     
                       r 
                       
                         1 
                          
                         
                           / 
                         
                          
                         m 
                       
                     
                   
                   〉 
                 
                 
                   r 
                   ≥ 
                   
                     
                       r 
                       min 
                     
                      
                     
                         
                     
                      
                     and 
                      
                     
                         
                     
                      
                     cell 
                      
                     
                         
                     
                      
                     not 
                      
                     
                         
                     
                      
                     in 
                      
                     
                         
                     
                      
                     tight 
                      
                     
                         
                     
                      
                     region 
                   
                 
               
               . 
             
           
         
       
     
         [0037]    Typically, method  15  is performed for a region of the borehole that is not larger than the resolution of the reference porosity log. The method  15  may be performed at multiple regions of interest, wherein the normalization factor and image point factors are calculated independently for each region of interest. 
         [0038]    At operation  26  of method  15 , the porosity is assigned based on the normalization factor and the image point factor. If the image point factor is smaller than the normalizing factor, equation (6) can be used to assign porosity to all cells which are not in the tight region (there the porosity is modeled as zero). If the image factor is not smaller than the normalizing factor, all cells of this resistivity are assigned to secondary porosity, and the minimum resistivity r min  for the matrix cells, the secondary porosity φ ν  and fractional volume V are updated accordingly, whereby the minimum matrix resistivity is set to the next lowest resistivity value. If the secondary porosity does not exceed a reasonable limit (e.g. 25%), referred to as the termination criteria in  FIG. 2 , operation  28 , the method proceeds to the next iteration, where the next lowest resistivity value is tested as to whether it satisfies inequality (8) which would qualify it for the lowest matrix resistivity. The example of 25% is not meant to be limiting; the termination criteria may be provided by the user based on any known or assumed properties of the geologic region of interest. An example is presence of very large vugs or caverns, and changes in caliper can be used to define a higher limit on secondary porosity in the termination criteria. The iterations are carried out until the inequality (8) is satisfied for r i  set to minimum resistivity in the matrix cells or the fractional volume ν exceeds the high limit. In the latter case, the region of tight rock can be adjusted to exclude a larger portion of the region of interest (e.g. in the case of delineating by cut-offs, the largest resistivity considered to this point to belong to matrix is now reassigned to the tight rock, and the iterative method of assessing secondary porosity can proceed starting from 0 again. If no resistivity is found to satisfy the inequality (8) with ν&lt;25%, for example, for the given tight rock volume, the iterative method can be carried out after successive adjustments to the tight rock volume until there are no more points with resistivity in the upper half of the resistivity image value range for the well. In such a case, the assessment is made from other log information as to whether this is e.g. a region of predominantly tight rock or of a high proportion of shale. Such an assessment then proceeds to assign the cell porosity to a constant for each of the two groups of cells, φ high   _   res.  and φ low   _   res.  (high resistivity and low resistivity group) according to whether there is clay or very tight rock in the region, observing that the porosity of the region has to match a reference porosity. 
         [0039]    The high- and low-resistivity group of cells can be adjusted in regions of tight rock or clay when a correction of image artifacts is performed. In one embodiment, such a correction is done from considerations of geometry. It uses screening of the points in resistivity image identified to be in regions of secondary porosity, to identify likely clay layers or noise. If resistivity of clay and tight rock is reasonably assumed to be constant and this constant is known, this correction would involve adjustment to image values themselves, prior to entering the image porosity calculation. In a region where there is only tight rock and secondary porosity present, secondary porosity for the region is obtained as 
         [0000]      φ v =φ−φ high   _   res.   ×N   high   /N  
 
         [0000]    where N high  is the number of cells with high resistivity and whereby φ high   _   res.  is set to the average value for tight rock in the well. 
         [0040]    Although method  15  of  FIG. 2  shows the calculation of the image point factor  24  occurring after the calculation of the normalization factor  22 , this is not meant to be limiting. These calculations may be done in any order or concurrently. 
         [0041]    It may also be desirable to calibrate the method  15 . This may be done by comparing the resistivity image log with core measurements. For example, in regions where there are multiple vugs of similar size displayed both in core and resistivity image, vug sizes in core can be used to deduce which points in image features belong to vugs and which may be conductive image artifacts. In one embodiment, such a calibration step corrects for artifacts by assigning a lower value to the local average porosity φ high  for secondary-porosity regions (e.g. vugs, mold size, etc.). 
         [0042]    Once it is deduced which cells belong to the secondary porosity, each matrix cell is assigned porosity according to equation (6). An example of the intermediate and final results of method  15  may be seen in  FIG. 5 . Here, the depth is represented in column  60  and the raw image is seen in column  62 . The calibrated image is in column  64 . The porosity image is in column  66 ; the light shades indicate areas with low porosity and dark shades indicate higher porosity. The porosity is also represented in columns  67  and  68 . Column  67  shows the porosity distribution and column  68  shows the average values for the secondary porosity (dark, short dashes), the total porosity (long dash, short dash), and the local matrix porosity (light gray, solid line). To obtain the total pore volume in matrix in the region of interest, the average local matrix porosity is multiplied by matrix volume. To obtain the total volume occupied by secondary porosity in the region of interest, secondary porosity is multiplied by the total volume of the region. 
         [0043]    When the method  15  terminates, it is possible to validate the results. This may be done, for example, by identifying regions of high conductivity from other well logs such as a gamma-ray log and/or a caliper log and elimination of those regions and all adjacent high-conductivity points from secondary porosity. 
         [0044]    A system  700  for performing the method  15  of  FIG. 2  is schematically illustrated in  FIG. 6 . The system includes a data source/storage device  70  which may include, among others, a data storage device or computer memory. The data source/storage device  70  may contain resistivity image log data. The data from data source/storage device  70  may be made available to a processor  72 , such as a programmable general purpose computer. The processor  72  is configured to execute computer modules that implement method  15 . These computer modules may include a normalization module  74  for calculating a normalization factor, an image point module  75  for calculating an image point factor, and a porosity module  76  for comparing the normalization factor and image point factor to determine where secondary porosity exists. These modules may be implemented more than once in an iterative manner. The system may include interface components such as user interface  79 . The user interface  79  may be used both to display raw data and processed data and to allow the user to select among options for implementing aspects of the method. By way of example and not limitation, the porosity distribution computed on the processor  72  may be displayed on the user interface  79 , stored on the data storage device or memory  70 , or both displayed and stored. 
         [0045]    While in the foregoing specification this invention has been described in relation to certain preferred embodiments thereof, and many details have been set forth for purpose of illustration, it will be apparent to those skilled in the art that the invention is susceptible to alteration and that certain other details described herein can vary considerably without departing from the basic principles of the invention. In addition, it should be appreciated that structural features or method steps shown or described in any one embodiment herein can be used in other embodiments as well.