- Email: [email protected]
Statistical analysis of jointed rock data
0148-9062/91 S3.00 + 0.00
Int. J. Rock Mech. Min. Sci. & Geomech. Absn. Vol. 28, No. 5, pp. 3X-382, 1991 Printed in Great Britain. All rights reserved
Copyright Q 1991 Pcrgamon Press pk
Statistical Analysis of Jointed Rock Data A. SHAPIRO? J. L. DELPORTS A methodology for the analysis of jointed rock data is described. The situation where a body of rock is divided into distinct blocks by families of roughly parallel fracture surfaces (joints) is considered. During the excavation of a tunnel through such rock, some of the blocks come into contact with the excavation face and may form so-called keyblocks. We address the problems of estimating the number of keyblocks per square metre of an excavation face, evaluating the distribution of keyblock volumes and estimating the probability that a suspected keyblock, part of which is observed as a closed polygon on an excavation face, may form a complete keyblock. Monte-Carlo methotis and simulation are used.
INTRODUCTION
STOCHASTIC MODELS
aim of this paper is to describe a methodology for statistical analysis of jointed rock data. We consider the situation where a body of rock is divided into distinct blocks by families (sets) of roughly parallel fracture surfaces (joints). During excavation of a tunnel through such rock, some of these blocks come into contact with the excavation face and may form so-called keyblocks. The problem of keyblock analysis and identification was discussed by Shi and Goodman [l, 21, Delport and Martin [3] and Vreede and Martin [4]. Here we concentrate on statistical aspects of keyblock analysis. In particular, the following problems will be addressed:
In this section, we discuss various statistical assump tions involved in the construction of an appropriate stochastic model. Borehole data were used to evaluate dip angles and dip directions of joints, and then to cluster them into distinct families (sets). The same data provided samples of fracture spacings corresponding to vertical distances between consecutive fractures from a particular family. Statistical distributions of these fracture spacings were analyzed. Let X,, . . . , X, be measurements of fracture spacings from a considered joint set. We view Xi, i = 1, . . . , n, as observations (a sample) of a random variable X representing the considered spacings. The empirical distributions of X have been compared with theoretical models suggested by the exponential and log-normal distributions. The exponential model is defined by the cumulative distribution function F(x; 2) = probability {X < x} depending on a single parameter rZ:
The
(i) estimating the number of keyblocks per square metre of an excavation face (tunnel wall); (ii) evaluating the probability distribution of keyblock volumes; and (iii) estimating the probability that a suspected keyblock, part of which is observed as a closed polygon on an excavation face, may form a complete keyblock. Problems (i) and (ii) were studied in McCullagh and Lang [5l where the exponential model for fracture spacings involving three fracture (joint) sets and infinite fracture trace lengths were employed. In this paper, we consider a more realistic model where the effects of the finiteness of joints are taken into account (cf. [6-S]). The complexity of problems thus obtained is such that a closed-form solution can hardly be expected. Therefore we resort to Monte-Carlo methods and simulation. Details of the statistical framework and associated algorithms are described in the following sections.
F(x;A)=
x >O.
In general the parameter A is unknown and should be estimated from the data. The standard maximum likelihood estimator of 1 is (see, e.g. [9]):
The log-normal distributions are characterized by two parameters p and tr with density function f(x;/f,a)=dF(x;/.4,a)/dx: f(x;p,o)=&exp
tuniversity of South Africa, P.O. Box 392, Pretoria 0001, South Africa. $Dcpartmcnt of Computational and Applied Mathematics, University of the Witwatcrsrand, Wits 2050, South Africa.
1 -exp(-rlx),
s(
- pg},
x
>o.
Notice that Y = In X transforms the log-normal variable X into a normal variable Y with mean p and standard 315
SHAPIRO and DELPORT:
316
STATISTICAL ANALYSIS OF JOINTED ROCK DATA
deviation u. Consequently the maximum likelihood estimators of p and 0 are: F = f ,f In Xi 1-I
Ouantlle
plot
for
first
data set
(2)
and
respectively. The expected (mean) value E(X) and the variance Var(X) of the log-normal variable X is related to the parameters p and LTas follows:
Var(X) = exp(2p + fr’} {exp(fr2) - l}.
(4)
Predicted
The expected value E(X) and the variance Var(X) can be estimated directly by:
Ouantae
2.4 t
plot for
spacing second
im)
data
set
(b)
and s2 =
& $I(Xi-X)‘, t=l
respectively. Then by solving equations (4) one obtains the following estimators of p and 6:
x
@ =ln
s2
i (
pf’
d =[ln{$+ L
\--
1;2
>
I ’
(5)
l}-/‘2.
Predicted
spacing (ml
/J
Estimators fi and d are consistent, although less efficient than the maximum likelihood estimators fi and 6. Theoretical distributions of the exponential and lognormal models have been plotted against empirical distributions. The quantile plots obtained for data sets of fracture spacings are given in Figs la,b and 2a,b. For the log-normal model the data have been transformed by the logarithm function and compared with the corresponding normal distribution. The straight line on the plots is a reference line with a slope of one. The more curved the empirical line, and the further away it is from the reference line, the stronger is the indication that the theoretical model does not fit the data. It can be seen that both models fit the data fairly well. The goodness of fit of the theoretical models has also been tested by the Kolmogorov-Smimov and Chisquare tests (see, e.g. [9]). In Table 1, n is the number of observations, d, is the sample value of the Kolmogorov-Smimov statistic D,, and the significance level is given by the probability Pr{D, > d,}. Notice that the parameters of exponential and log-normal distributions were estimated from the data. Therefore the results of Kolmogorov-Smimov test should not be taken literally (Table 2).
Fig. I
It appears that both models fit the data reasonably well with a slight trend in favour of the exponential model. Recall that the significance level is a number between 0 and l-the closer to 1, the better the fit. Moreover, the exponential model has a number of important theoretical properties which make it more suitable for a simulation implementation. This will be discussed further in the next section. The intersection of joints with an excavation surface results in line segments visible on the corresponding excavation face. The statistical distribution of the lengths of these line segments (trace lengths) is the second major factor in construction of the model. Unfortunately no data sets were available to us and it was necessary to rely on studies found in the literature. In this respect we follow Rouleau and Gale [8] where it is suggested that the log-normal model provides a good fit for the data. Given a sample 2, , . . . , ZN of trace lengths (supposedly coming from the same set of observations as spacings), the required parameters ~1and 0 can be estimated by employing the maximum likelihood statistics fi and 6, given in (2) and (3), respectively, where X, is replaced by Zi and n by N. If the sample mean 2 and
Ouantilr
plot
for
data
first
Table 1. Results of the Kolmogorov-Smimov test for the exponential :and log-normal models
ret
Log-normal
Exponential
(a)
l
Line
317
STATISTICAL ANALYSIS OF JOINTED ROCK DATA
SHAPIRO and DELPORT:
Set 1 Set 2
W=SZ
+ Lop-normal
n
d.
29 31
0.083 0.099
Significance level 1
1
n
dn
29 31
0.084 0.110
Significance level 1
1
will generate spacings between parallel line segments in accordance with a given exponential distribution; and (ii) lengths of line segments from the same family have a log-normal distribution with given parameters p and Q.
Predicted
Ouantile
logarithm
plot
for second
of
In order to proceed, the following observation must be made. Consider a sequence xl, . . . . x, of random numbers uniformly and independently distributed on the interval [0,11.Let xc,) d x(r) < . . . G xc,) be this sequence rearranged in increasing order (called order statistics), i.e. xu, is the smallest of the numbers xi, i = 1, . . . , n, etc. Then for all i the length xci+,) -xc,) of ith interval has the same distribution as the smallest variable xc,) (e.g. Kendall and Moran [lo, $2.81). The cumulative distribution function F(x) of xcI), and hence of every generated interval, is:
spacing
data
ret
2
(b) t
F(x)=!-
l-5 “9 OGx,
I -1.9
I -0.9 Predicted
I
I
0.1
1.1
logarithm
of
I 2 .l
spacing
Fig. 2
variance S2 are provided, then p and cr can be estimated by ji and d, given in (5) and (6), respectively. MODEL CONSTRUCTION
in this section we outline a procedure for simulation of an excavation face. Families of parallel line segments, representing traces of joints, should be constructed in such a way that: (i) an imaginary borehole, drawn at random and perpendicular to the base of a rectangular picture,
This shows that lengths of the generated intervals are distributed approximately exponentially with the parameter 3, = n/l. (This result can be derived rigorously. It can be shown that for large n the asymptotic distribution of “xc, is exponential with the corresponding cumulative distribution function, 1 - exp( - x/l), e.g. [1 1, Section 2.3. I].) Now consider a rectangle with the base of length L and height H. Let us generate random numbers , y,., , uniformly and independently XI, - * *, xNandy,,... distributed on the intervals [0, L] and [0, H], respectively, and form N points in the rectangle with coordinates (xi, yi), i = 1, . . . , N. Furthermore, generate random numbers z,, . . . , zN independently distributed according to a log-normal distribution with given parameters p and 6. [This can be done by generating random numbers wI, . . . , wN, distributed normally N@, a), and taking Zi = exp(wi).] For each generated point (x,, n) draw the line segment of length zi centred at (x,, J+) and at a given angle $ to the base of the
Table 2. Results of the &i-square test for the exponential and log-normal models Exponential
Set I Set 2 UlMS
2dls-c
Log-normal
n
x2
d.f.
Sigaifieanee level
29 31
0.854 2.675
2 3
0.653 0.444
z exp(-nx/l).
n
x2
d.f.
29 31
2.099 2.743
2
1
Significance level 0.147 0.254
SHAPIRO and DELPORT:
378
STATISTICAL ANALYSIS OF JOINTED ROCK DATA x
length of the segment, starting the rectangle
from this
length zi of the considered ated. An example of the picture obtained
gener-
I
I
I
I b
second coordinate
w
Fig. 4
inside the rectangle. p of such intersection. Let Z be a log-normal random variable, W = In Z be the corresponding normal variable iV(p, a) and X be a random variable distributed uniformly on the interval [0, L]. The variables Z and X are assumed to be independent. Density functionsf&w) and fX(x) of W and X are:
respectively. Then: p = Pr{Zlcos $ I 2 X} = Pr{ew]cos $1 2 X},
where 0 6 ((I < x. This probability integral: P=
,(w)=*exp{-@&$}
is given by the
fx(x)_fw(w) dx dw, ss A
where A is the region
and
A = {(x, w): ew]cos + I 2 x, 0 < x < L}.
osx
fx(x)=L-‘,
GL,
Let b be the number such that lcos 1(1]exp(b) = L, i.e. b =lnL
-ln]costj].
Then 1 P=La*
b
ew1cos $1
s --ccs 0
xexp{ -w}
=-
1 Lq/G
Length
of borehole
13.00m
Length
of exposure
6.00m
Height
01 exposure
6.23m
Height
01 tunnel
2.00m
Dip of excavation DIP direction Fault
: 3 4
No on
Dip
Bearmg
No in borehole
65”
145”
36
55
f;:
::o 346’
: 2
z 5
70”
0
exposure
Fig. 3
a(x)
110.00’
face
Tracelength mean
STD
3.0m
l.Om
lcos$]exp{w
-&$ldw
s --cc
where
90.00’
face
of exCavatiOn
b
dxdw +Pr(W>b}
Probability
angle 0.3154 0.4809 0.2541 0.2010
=-&
lI_exp(
-G)d,
Base 129” 23’
is the cumulative distribution normal variable. Now (w - p)*
w---ZT=-
function of a standard
[w - (p + 2a2
a2)12 + ~ +;
SHAPIRO and DELPORI:
379
STATISTICAL ANALYSIS OF JOINTED ROCK DATA
given family, in this portion of the borehole. Therefore in the case of a vertical excavation face:
and hence expIw - (w2-;)‘)
N =p-‘m(H/i). - kw -;C:e*))2)exp(/l
= exp
For a non-vertical excavation face (wall) the following correction must be made. Let II = (n,, nY,n,) be a unit normal to joints of the considered family, and u = (u,, u,,, u,) be a unit vector in the direction corresponding to the straight lines perpendicular to the base of the generated picture. For example, if n,, is a normal to the considered excavation face and a is a vector in the direction of the tunnel axis, then u can be taken perpendicular to a, and a, that is:
+$).
I Then
x exp
=exp
_ lw - (p + a2)12 2az
02 p +y
(10)
For a vertical excavation face, u is usually equal to the coordinate vector k = (O,O, 1). When the roof is considered, u is given by a unit normal to any one of the tunnel walls. Now the average number of intersections must be corrected by the factor In-kl-’ In*ul and hence:
( 1
Pr(U<6},
where U is normal N(p + rr2, a). Since
N = p-‘m(H/f)ln,l-‘(naui.
Notice that the probability p in (11) also depends on the choice of the excavation face through the corresponding angle $. Finally, we note that if the dip angle a and dip direction (bearing) /3 of an excavation face are given, then the coordinates of the corresponding unit normal are:
we finally obtain:
+1-Q
3 (
d
. >
(7)
The number p in (7) gives the probability that a generated line segment intersects the considered straight line only once. Because of reflections it may happen that, if the trace length is large enough, there is more than one intersection corresponding to the same line segment. To avoid this problem we restrict the generated trace lengths zi by taking: z:=min{zi,&}. That is, if the length z,Jcos I++ 1 of the line segment projection onto the rectangle base is greater than L, then we restrict it to L. We know that the expected number of intersections of the considered straight line with the generated line segments is Np, whereas it should be IZH. This suggests that the number N of generated line segments should be taken such that Np = AH, i.e. N =p-‘AH.
(11)
(9)
Of course, N must be an integer. Therefore the integer part of the right-hand side in (9) is taken. For a vertical excavation face, the considered straight line corresponds to a vertical borehole. In this case the parameter 1 is estimated by m/l, where I is the length of a considered portion of the borehole and m is the number of appearances of joints, corresponding to a
n,= sina sin/I, nY
=
sina cos/3,
n,=cosa. A portion of a vertical excavation face is simulated in Fig. 3. Two thick parallel lines in the centre of the picture mark an imaginary tunnel wall. There are four families of joints, traces of which form angles JI with the base of the picture. The “length of borehole” stands for the effective domain 1 of the borehole where measurements have been taken. The corresponding number m of appearances of joints is given under the heading “No. in borehole”. The number N of generated line segments is given under the heading “No. on exposure”. Under the heading “probability” are given probabilities p calculated according to the formula (7). The “mean” and “STD” denote the sample mean Z and sample standard deviation S of trace lengths. The corresponding parameters p and Q are estimated by ji and 5, given in (5) and (6), respectively. ESTIMATION
OF KEYBLOCK
PARAMETERS
The most frequent and therefore important case of keyblock occurrence involves three families of joints and an excavation face. This is the simplest and most manageable non-trivial situation to which we limit our attention in this section. Of course, the developed pro-
380
SHAPlRO and DELPORT:
STATlST7CAL ANALYSIS OF JOINTED ROCK DATA
cedure can be repeated for all possible triples of the considered families. Some of the generated line segments corresponding to the triple of chosen families form closed triangles. All such triangles are identified and subsequently tested by the keyblock procedure [3]. The retained triangles form bases of potential keyblocks. In Fig. 5 the simulation picture of Fig. 3 is analyzed for keyblocks associated with families 1, 2 and 3. The shaded triangles in Fig. 5 mark bases of possible or suspected keyblocks. It is of particular interest and practical importance to evaluate the probability that the corresponding joints extend deeply enough into the rock to form a complete keyblock. This problem will be discussed later. The simulation procedure is repeated several times in the computer memory until N = 300 keyblocks are created and counted. This provides an estimate of the number of suspected keyblocks per square metre on the excavation face. Notice that there are a considerable number of overlapping keyblocks (shaded differently) on the picture. Each one is counted separately in the estimation of “keyblocks per square metre”. The estimated number of keyblocks per square metre depends on the input parameters. In particular, the “length of exposure” and “height of exposure” represent the size of a portion of the excavation face under consideration. Since only the suspected keyblocks
which lie ensirefy wlthin this portion of the excavation face are counted, for larger values of these two parameters one would expect bigger values of the estimated number of keyblocks per square metre. The described estimation procedure was repeated 10 times for 10 x 10 and 15 x 15 (in m) pictures, with the remaining parameters given in Fig. 5. The obtained estimates of “keyblocks per square metre” were 1.02,0.88,0.73,0.82, 0.76,0.94,0.72, 0.76, 0.62. 0.92 (av. 0.82) and 0.73,0.88, 0.94, 0.91, 0.84, 0.71, 0.77, 0.69, 0.77, 1.12 (av. 0.84), respectively. Naturally, limitations of the tunnel size (height) will not allow for keyblocks that are “too big” to appear on the tunnel wall. To calculate the relevant estimates one can apply the described procedure when the “height of exposure” parameter is equal to the “height of tunnel”. The estimates obtained (for the case given in Fig. 5) when “height of exposure” and “length of exposure” are equal to 2 and 10 m, respectively, are 0.63, 0.60, 0.70, 0.58, 0.45, 0.52, 0.63, 0.57, 0.52, 0.58 (av. 0.58). When the “mean” and “STD” parameters were changed from 3 and 1 to 6 and 2, respectively, the procedure produced the following estimates for keyblocks per square metre (with “height of exposure” = 2 and “length of exposure” = 10): 0.82, 0.78, 0.68, 0.91, 0.84, 0.71, 0.78, 0.67, 0.67, 0.63 (av. 0.75). Similarly for “mean” = 1.5 and “STD” = 0.5 the estimates obtained were 0.49, 0.36, 0.38, 0.36, 0.36, 0.34, 0.41, 0.42, 0.37, 0.41, 0.38 (av. 0.39). Now in order to evaluate the distribution of keyblock volumes we proceed as follows: Letn,=(ni,,niZ,ni3), i = 1, 2, 3, beanormal tojointsof the ith family and n, = (no,, n,, no3) be a normal to the excavation face. Consider the pyramid bounded by the planes: njfX
+nfzY
+nfjZ
=O,
i = 1,2,3, and the plane
n,,x+n,,y+n,,z=l. This pyramid has a vertex at 0 = (0,0,O)and three other vertices v, = (xi, yi, z,), i = 1,2, 3, which can be calculated by solving three systems of linear equations. We can also calculate the angles yi of the pyramid adjacent to the vertex 0: 13.00m
Length of borehola Length
ot
Height
of exposure
6.23m
Height
of tunnel
2.00m
Dip of excavation Dip direction Fault
Dip
(12)
6.00m
exposure
Searing
90.00”
face
of excavation
No in borehole --36 6 :
face
No on exooswe 55 6 65
Fig. 5
110.00”
Tracelength mean STD 3.0m 1.0m 3.0m l.Om 3.0m l.Om
Probability 0.3154 0.4609 0.2010 0.2541
Base anale 129’ 23” :to 0
i=1,2,3. All considered keyblocks are similar to each other and to the pyramid constructed above. Therefore y, represents the angles of these pyramids adjacent to the vertex inside the rock. Moreover, let S be the area of the pyramid base visible on the excavation face. Then the volume V of this pyramid is: V = KS3’2,
(13)
SHAPIRO and DELPORTz Cumulative distribution
I
0
0.04
381
of keyblock volumrs
I
I
0.12
0.16
I
0.06
STATISTICAL ANALYSIS OF JOINTED ROCK DATA
I
0.20
I
0.24
Volume Fig. 6
where the coefficient K is the same for all considered keyblock pyramids and can be calculated as follows: Consider the constructed pyramid with vertices v, , v2, v3 and 0. The volume V* of this pyramid is: XI
det i
YI
ZI
x2
Y2
z2
x3
Y3
z3
1
and the height h* is 1
h* =.
_ JG,
_ + nt
_* + nb
(The height h * is equal to one if ]ln, 1)= 1.) Then: I/* = (fh*)“z(V*)‘/’
(14)
Formula (13) allows quick calculation of keyblock volumes. Volumes of the identified keyblocks, simulated by the procedure, were calculated and stored in the computer memory. The cumulative distribution of these volumes was then calculated. An example of such distribution (for the tunnel wall of the case described in Fig. 5) is shown in Fig. 6. The cumulative distribution of Fig. 6 indicates, for instance, that about 90% of the keyblocks have a volume less than 0.06 m3. In the remainder of this section the problem of estimating the probability that a suspected keyblock may form a complete keyblock is discussed. Consider a line segment AB which forms side QR of the triangle AQPR representing the base of a suspected keyblock. The line segment AB represents the trace of intersection of a joint with the excavation face. This joint forms
Fig. 8
a face QRC of the keyblock pyramid with the vertex C inside the rock. The shape of joint is idealized to form a circle and it is assumed that the plane of the excavation face intersects them at “random”. We have to calculate the probability pi that the triangle AQRC constructed on the interval QR of the “random” chord AB falls entirely within the circle. Of course, this probability depends on the distribution of the “random” chord AB and distances AQ = bi, QR = c, and RB = di. We assume that the intersection of the excavation plane and the joint is arbitrary, and that therefore the intersection D of AB with the perpendicular diameter is uniformly distributed on that diameter. Then pi is given by the proportion of the diameter cut off by the chord AB when the point C lies on the circumference of the circle. After certain calculations we obtain: pi = [sin(8,/2)12, where Ji is the angle L ACB which can be calculated as follows. Consider the angles:
A 0
Fig. I
Fig. 9
SHAIIIP.0 and DELPORT:
?82
STATJSTJCAL ANALYSIS OF JOINTED ROCK DATA
the same family. However, by definition, such a quadrangle cannot be a base of a keyblock (see [3]). Extension of the estimation procedure to keyblocks formed by families of four or more joints, is a matter for further development. Notice that once the model simulation, described in the third section, is completed it basically becomes the question of computer time to count ail simulated keyblocks formed by families of, say, four joints and to calculate the corresponding estimates.
and /I& =
L
CRQ =
7T -
a, --
3’, .
Then
where sin ii = dJ(sin a,/sin yi)‘cf + df - 2cidi(sin a,/sin y,jcos(a, + yi)]-“2 Sin(a, + ri)
Acknowledgements-This paper is based on a project supported by the Centre for Advanced Computing and Decision Support, CSIR, Pretoria and follows CSIR Technical Report PKOMP 89/2. The authors wish to express their thanks to H. W. Ittman for careful management of the project, to F. A. Vreede for helpful suggestions and to the referee for constructive comments.
and sin vi = b,[(sin /$/sin yi)2cf + bf - 2cibi(sin /Ii/Sin r,)COS(/3i+ yi)]-“’ Sin(fii + ri). If we assume that joints of different families are independently distributed, then the probability P that a suspected keyblock forms a complete keyblock is p,p2pj. It follows that under the above model: P = [sin(s,/2)sin(82/2)sin(6,/2)]2.
Notice that the angles ai and hence the probability depend on the distances bi, Ci, 4, i = 1,2,3. CONCLUDING
(15) P
REMARKS
It should be understood that the procedure developed here gives only an approximation of reality. Large errors can occur due to inaccurate estimates of the parameters p and r~ of the log-normal distribution of trace lengths, because for instance, line segments of an identified family are not exactly parallel. It is recommended that the program should be run several times with slightly different values of the input parameters, in order to get some idea about variability of the output estimators. It is emphasized again that all suspected keyblocks visible on the simulated excavation face should be counted even if they overlap each other (which is often the case). Notice that it is quite possible for line segments from three families under consideration to form a closed quadrangle with two sides formed by parallel lines from
Accepted for publication 25 January 1991.
REFERENCES 1. Shi G. H. and Goodman R. E. Underground support design using block theory to determine keyblock bolting requirements. Proc. SANGORM
Symp. on Rock Mech.
in the Design of Tunnels,
Pretoria, pp. 81-105 (1963). 2. Shi G. H. and Goodman R. E. Keyblock bolting. Proc. Inr. Symp. on Rock Bolting, Abisko, pp. 143-167 (1983). 3. Delport J. L. and Martin D. H. A multiplier method for identifying keyblocks in excavations through jointed rock. SIAM J. Aig. Disc. Meth. 7, 321-330 (1986). 4. Vreede F. A. and Martin D. H. Keyblock characterization of jointed rock. SANGORM Symp. on Rock Mass Characterization (1985).
5. McCullagh P. and Lang P. Stochastic models for rock instability in tunnels. J. R. Statist. Sot. B. 46, 344-352 (1984). 6. Priest S. D. and Hudson J. A. Estimation of discontinuity spacing and trace length using scanline surveys. ht. J. Rock. Mech. Mm. Sci. & Geomech. Abstr. 18, 183-197 (1981).
7. Pahl P. J. Estimating the mean length of discontinuity traces. Inf. J. Rock Mech. Min. Sci. & Geomech. Abstr. X3, 221-228 (1981). 8. Rouleau A. and Gale J. E. Statistical characterization of the fracture system in the Stripa granite. Sweden. Int. J. Rock Mech. Min. Sci. & Geomech. Abstr. 22, 353-367 (1985). 9. Bury K. V. Statistical Models in Applied Sciences. Wiley, New York (1975). 10. Kendall M. G. and Moran P. A. P. Geometrical Probability. Griffin’s Statistical Monographs, London (1963). 11. Galambos J. The Asymptotic Theory of Extreme Order Statistics. Wiley, New York (1978).
Int. J. Rock Mech. Min. Sci. & Geomech. Absn. Vol. 28, No. 5, pp. 3X-382, 1991 Printed in Great Britain. All rights reserved
Copyright Q 1991 Pcrgamon Press pk
Statistical Analysis of Jointed Rock Data A. SHAPIRO? J. L. DELPORTS A methodology for the analysis of jointed rock data is described. The situation where a body of rock is divided into distinct blocks by families of roughly parallel fracture surfaces (joints) is considered. During the excavation of a tunnel through such rock, some of the blocks come into contact with the excavation face and may form so-called keyblocks. We address the problems of estimating the number of keyblocks per square metre of an excavation face, evaluating the distribution of keyblock volumes and estimating the probability that a suspected keyblock, part of which is observed as a closed polygon on an excavation face, may form a complete keyblock. Monte-Carlo methotis and simulation are used.
INTRODUCTION
STOCHASTIC MODELS
aim of this paper is to describe a methodology for statistical analysis of jointed rock data. We consider the situation where a body of rock is divided into distinct blocks by families (sets) of roughly parallel fracture surfaces (joints). During excavation of a tunnel through such rock, some of these blocks come into contact with the excavation face and may form so-called keyblocks. The problem of keyblock analysis and identification was discussed by Shi and Goodman [l, 21, Delport and Martin [3] and Vreede and Martin [4]. Here we concentrate on statistical aspects of keyblock analysis. In particular, the following problems will be addressed:
In this section, we discuss various statistical assump tions involved in the construction of an appropriate stochastic model. Borehole data were used to evaluate dip angles and dip directions of joints, and then to cluster them into distinct families (sets). The same data provided samples of fracture spacings corresponding to vertical distances between consecutive fractures from a particular family. Statistical distributions of these fracture spacings were analyzed. Let X,, . . . , X, be measurements of fracture spacings from a considered joint set. We view Xi, i = 1, . . . , n, as observations (a sample) of a random variable X representing the considered spacings. The empirical distributions of X have been compared with theoretical models suggested by the exponential and log-normal distributions. The exponential model is defined by the cumulative distribution function F(x; 2) = probability {X < x} depending on a single parameter rZ:
The
(i) estimating the number of keyblocks per square metre of an excavation face (tunnel wall); (ii) evaluating the probability distribution of keyblock volumes; and (iii) estimating the probability that a suspected keyblock, part of which is observed as a closed polygon on an excavation face, may form a complete keyblock. Problems (i) and (ii) were studied in McCullagh and Lang [5l where the exponential model for fracture spacings involving three fracture (joint) sets and infinite fracture trace lengths were employed. In this paper, we consider a more realistic model where the effects of the finiteness of joints are taken into account (cf. [6-S]). The complexity of problems thus obtained is such that a closed-form solution can hardly be expected. Therefore we resort to Monte-Carlo methods and simulation. Details of the statistical framework and associated algorithms are described in the following sections.
F(x;A)=
x >O.
In general the parameter A is unknown and should be estimated from the data. The standard maximum likelihood estimator of 1 is (see, e.g. [9]):
The log-normal distributions are characterized by two parameters p and tr with density function f(x;/f,a)=dF(x;/.4,a)/dx: f(x;p,o)=&exp
tuniversity of South Africa, P.O. Box 392, Pretoria 0001, South Africa. $Dcpartmcnt of Computational and Applied Mathematics, University of the Witwatcrsrand, Wits 2050, South Africa.
1 -exp(-rlx),
s(
- pg},
x
>o.
Notice that Y = In X transforms the log-normal variable X into a normal variable Y with mean p and standard 315
SHAPIRO and DELPORT:
316
STATISTICAL ANALYSIS OF JOINTED ROCK DATA
deviation u. Consequently the maximum likelihood estimators of p and 0 are: F = f ,f In Xi 1-I
Ouantlle
plot
for
first
data set
(2)
and
respectively. The expected (mean) value E(X) and the variance Var(X) of the log-normal variable X is related to the parameters p and LTas follows:
Var(X) = exp(2p + fr’} {exp(fr2) - l}.
(4)
Predicted
The expected value E(X) and the variance Var(X) can be estimated directly by:
Ouantae
2.4 t
plot for
spacing second
im)
data
set
(b)
and s2 =
& $I(Xi-X)‘, t=l
respectively. Then by solving equations (4) one obtains the following estimators of p and 6:
x
@ =ln
s2
i (
pf’
d =[ln{$+ L
\--
1;2
>
I ’
(5)
l}-/‘2.
Predicted
spacing (ml
/J
Estimators fi and d are consistent, although less efficient than the maximum likelihood estimators fi and 6. Theoretical distributions of the exponential and lognormal models have been plotted against empirical distributions. The quantile plots obtained for data sets of fracture spacings are given in Figs la,b and 2a,b. For the log-normal model the data have been transformed by the logarithm function and compared with the corresponding normal distribution. The straight line on the plots is a reference line with a slope of one. The more curved the empirical line, and the further away it is from the reference line, the stronger is the indication that the theoretical model does not fit the data. It can be seen that both models fit the data fairly well. The goodness of fit of the theoretical models has also been tested by the Kolmogorov-Smimov and Chisquare tests (see, e.g. [9]). In Table 1, n is the number of observations, d, is the sample value of the Kolmogorov-Smimov statistic D,, and the significance level is given by the probability Pr{D, > d,}. Notice that the parameters of exponential and log-normal distributions were estimated from the data. Therefore the results of Kolmogorov-Smimov test should not be taken literally (Table 2).
Fig. I
It appears that both models fit the data reasonably well with a slight trend in favour of the exponential model. Recall that the significance level is a number between 0 and l-the closer to 1, the better the fit. Moreover, the exponential model has a number of important theoretical properties which make it more suitable for a simulation implementation. This will be discussed further in the next section. The intersection of joints with an excavation surface results in line segments visible on the corresponding excavation face. The statistical distribution of the lengths of these line segments (trace lengths) is the second major factor in construction of the model. Unfortunately no data sets were available to us and it was necessary to rely on studies found in the literature. In this respect we follow Rouleau and Gale [8] where it is suggested that the log-normal model provides a good fit for the data. Given a sample 2, , . . . , ZN of trace lengths (supposedly coming from the same set of observations as spacings), the required parameters ~1and 0 can be estimated by employing the maximum likelihood statistics fi and 6, given in (2) and (3), respectively, where X, is replaced by Zi and n by N. If the sample mean 2 and
Ouantilr
plot
for
data
first
Table 1. Results of the Kolmogorov-Smimov test for the exponential :and log-normal models
ret
Log-normal
Exponential
(a)
l
Line
317
STATISTICAL ANALYSIS OF JOINTED ROCK DATA
SHAPIRO and DELPORT:
Set 1 Set 2
W=SZ
+ Lop-normal
n
d.
29 31
0.083 0.099
Significance level 1
1
n
dn
29 31
0.084 0.110
Significance level 1
1
will generate spacings between parallel line segments in accordance with a given exponential distribution; and (ii) lengths of line segments from the same family have a log-normal distribution with given parameters p and Q.
Predicted
Ouantile
logarithm
plot
for second
of
In order to proceed, the following observation must be made. Consider a sequence xl, . . . . x, of random numbers uniformly and independently distributed on the interval [0,11.Let xc,) d x(r) < . . . G xc,) be this sequence rearranged in increasing order (called order statistics), i.e. xu, is the smallest of the numbers xi, i = 1, . . . , n, etc. Then for all i the length xci+,) -xc,) of ith interval has the same distribution as the smallest variable xc,) (e.g. Kendall and Moran [lo, $2.81). The cumulative distribution function F(x) of xcI), and hence of every generated interval, is:
spacing
data
ret
2
(b) t
F(x)=!-
l-5 “9 OGx,
I -1.9
I -0.9 Predicted
I
I
0.1
1.1
logarithm
of
I 2 .l
spacing
Fig. 2
variance S2 are provided, then p and cr can be estimated by ji and d, given in (5) and (6), respectively. MODEL CONSTRUCTION
in this section we outline a procedure for simulation of an excavation face. Families of parallel line segments, representing traces of joints, should be constructed in such a way that: (i) an imaginary borehole, drawn at random and perpendicular to the base of a rectangular picture,
This shows that lengths of the generated intervals are distributed approximately exponentially with the parameter 3, = n/l. (This result can be derived rigorously. It can be shown that for large n the asymptotic distribution of “xc, is exponential with the corresponding cumulative distribution function, 1 - exp( - x/l), e.g. [1 1, Section 2.3. I].) Now consider a rectangle with the base of length L and height H. Let us generate random numbers , y,., , uniformly and independently XI, - * *, xNandy,,... distributed on the intervals [0, L] and [0, H], respectively, and form N points in the rectangle with coordinates (xi, yi), i = 1, . . . , N. Furthermore, generate random numbers z,, . . . , zN independently distributed according to a log-normal distribution with given parameters p and 6. [This can be done by generating random numbers wI, . . . , wN, distributed normally N@, a), and taking Zi = exp(wi).] For each generated point (x,, n) draw the line segment of length zi centred at (x,, J+) and at a given angle $ to the base of the
Table 2. Results of the &i-square test for the exponential and log-normal models Exponential
Set I Set 2 UlMS
2dls-c
Log-normal
n
x2
d.f.
Sigaifieanee level
29 31
0.854 2.675
2 3
0.653 0.444
z exp(-nx/l).
n
x2
d.f.
29 31
2.099 2.743
2
1
Significance level 0.147 0.254
SHAPIRO and DELPORT:
378
STATISTICAL ANALYSIS OF JOINTED ROCK DATA x
length of the segment, starting the rectangle
from this
length zi of the considered ated. An example of the picture obtained
gener-
I
I
I
I b
second coordinate
w
Fig. 4
inside the rectangle. p of such intersection. Let Z be a log-normal random variable, W = In Z be the corresponding normal variable iV(p, a) and X be a random variable distributed uniformly on the interval [0, L]. The variables Z and X are assumed to be independent. Density functionsf&w) and fX(x) of W and X are:
respectively. Then: p = Pr{Zlcos $ I 2 X} = Pr{ew]cos $1 2 X},
where 0 6 ((I < x. This probability integral: P=
,(w)=*exp{-@&$}
is given by the
fx(x)_fw(w) dx dw, ss A
where A is the region
and
A = {(x, w): ew]cos + I 2 x, 0 < x < L}.
osx
fx(x)=L-‘,
GL,
Let b be the number such that lcos 1(1]exp(b) = L, i.e. b =lnL
-ln]costj].
Then 1 P=La*
b
ew1cos $1
s --ccs 0
xexp{ -w}
=-
1 Lq/G
Length
of borehole
13.00m
Length
of exposure
6.00m
Height
01 exposure
6.23m
Height
01 tunnel
2.00m
Dip of excavation DIP direction Fault
: 3 4
No on
Dip
Bearmg
No in borehole
65”
145”
36
55
f;:
::o 346’
: 2
z 5
70”
0
exposure
Fig. 3
a(x)
110.00’
face
Tracelength mean
STD
3.0m
l.Om
lcos$]exp{w
-&$ldw
s --cc
where
90.00’
face
of exCavatiOn
b
dxdw +Pr(W>b}
Probability
angle 0.3154 0.4809 0.2541 0.2010
=-&
lI_exp(
-G)d,
Base 129” 23’
is the cumulative distribution normal variable. Now (w - p)*
w---ZT=-
function of a standard
[w - (p + 2a2
a2)12 + ~ +;
SHAPIRO and DELPORI:
379
STATISTICAL ANALYSIS OF JOINTED ROCK DATA
given family, in this portion of the borehole. Therefore in the case of a vertical excavation face:
and hence expIw - (w2-;)‘)
N =p-‘m(H/i). - kw -;C:e*))2)exp(/l
= exp
For a non-vertical excavation face (wall) the following correction must be made. Let II = (n,, nY,n,) be a unit normal to joints of the considered family, and u = (u,, u,,, u,) be a unit vector in the direction corresponding to the straight lines perpendicular to the base of the generated picture. For example, if n,, is a normal to the considered excavation face and a is a vector in the direction of the tunnel axis, then u can be taken perpendicular to a, and a, that is:
+$).
I Then
x exp
=exp
_ lw - (p + a2)12 2az
02 p +y
(10)
For a vertical excavation face, u is usually equal to the coordinate vector k = (O,O, 1). When the roof is considered, u is given by a unit normal to any one of the tunnel walls. Now the average number of intersections must be corrected by the factor In-kl-’ In*ul and hence:
( 1
Pr(U<6},
where U is normal N(p + rr2, a). Since
N = p-‘m(H/f)ln,l-‘(naui.
Notice that the probability p in (11) also depends on the choice of the excavation face through the corresponding angle $. Finally, we note that if the dip angle a and dip direction (bearing) /3 of an excavation face are given, then the coordinates of the corresponding unit normal are:
we finally obtain:
+1-Q
3 (
d
. >
(7)
The number p in (7) gives the probability that a generated line segment intersects the considered straight line only once. Because of reflections it may happen that, if the trace length is large enough, there is more than one intersection corresponding to the same line segment. To avoid this problem we restrict the generated trace lengths zi by taking: z:=min{zi,&}. That is, if the length z,Jcos I++ 1 of the line segment projection onto the rectangle base is greater than L, then we restrict it to L. We know that the expected number of intersections of the considered straight line with the generated line segments is Np, whereas it should be IZH. This suggests that the number N of generated line segments should be taken such that Np = AH, i.e. N =p-‘AH.
(11)
(9)
Of course, N must be an integer. Therefore the integer part of the right-hand side in (9) is taken. For a vertical excavation face, the considered straight line corresponds to a vertical borehole. In this case the parameter 1 is estimated by m/l, where I is the length of a considered portion of the borehole and m is the number of appearances of joints, corresponding to a
n,= sina sin/I, nY
=
sina cos/3,
n,=cosa. A portion of a vertical excavation face is simulated in Fig. 3. Two thick parallel lines in the centre of the picture mark an imaginary tunnel wall. There are four families of joints, traces of which form angles JI with the base of the picture. The “length of borehole” stands for the effective domain 1 of the borehole where measurements have been taken. The corresponding number m of appearances of joints is given under the heading “No. in borehole”. The number N of generated line segments is given under the heading “No. on exposure”. Under the heading “probability” are given probabilities p calculated according to the formula (7). The “mean” and “STD” denote the sample mean Z and sample standard deviation S of trace lengths. The corresponding parameters p and Q are estimated by ji and 5, given in (5) and (6), respectively. ESTIMATION
OF KEYBLOCK
PARAMETERS
The most frequent and therefore important case of keyblock occurrence involves three families of joints and an excavation face. This is the simplest and most manageable non-trivial situation to which we limit our attention in this section. Of course, the developed pro-
380
SHAPlRO and DELPORT:
STATlST7CAL ANALYSIS OF JOINTED ROCK DATA
cedure can be repeated for all possible triples of the considered families. Some of the generated line segments corresponding to the triple of chosen families form closed triangles. All such triangles are identified and subsequently tested by the keyblock procedure [3]. The retained triangles form bases of potential keyblocks. In Fig. 5 the simulation picture of Fig. 3 is analyzed for keyblocks associated with families 1, 2 and 3. The shaded triangles in Fig. 5 mark bases of possible or suspected keyblocks. It is of particular interest and practical importance to evaluate the probability that the corresponding joints extend deeply enough into the rock to form a complete keyblock. This problem will be discussed later. The simulation procedure is repeated several times in the computer memory until N = 300 keyblocks are created and counted. This provides an estimate of the number of suspected keyblocks per square metre on the excavation face. Notice that there are a considerable number of overlapping keyblocks (shaded differently) on the picture. Each one is counted separately in the estimation of “keyblocks per square metre”. The estimated number of keyblocks per square metre depends on the input parameters. In particular, the “length of exposure” and “height of exposure” represent the size of a portion of the excavation face under consideration. Since only the suspected keyblocks
which lie ensirefy wlthin this portion of the excavation face are counted, for larger values of these two parameters one would expect bigger values of the estimated number of keyblocks per square metre. The described estimation procedure was repeated 10 times for 10 x 10 and 15 x 15 (in m) pictures, with the remaining parameters given in Fig. 5. The obtained estimates of “keyblocks per square metre” were 1.02,0.88,0.73,0.82, 0.76,0.94,0.72, 0.76, 0.62. 0.92 (av. 0.82) and 0.73,0.88, 0.94, 0.91, 0.84, 0.71, 0.77, 0.69, 0.77, 1.12 (av. 0.84), respectively. Naturally, limitations of the tunnel size (height) will not allow for keyblocks that are “too big” to appear on the tunnel wall. To calculate the relevant estimates one can apply the described procedure when the “height of exposure” parameter is equal to the “height of tunnel”. The estimates obtained (for the case given in Fig. 5) when “height of exposure” and “length of exposure” are equal to 2 and 10 m, respectively, are 0.63, 0.60, 0.70, 0.58, 0.45, 0.52, 0.63, 0.57, 0.52, 0.58 (av. 0.58). When the “mean” and “STD” parameters were changed from 3 and 1 to 6 and 2, respectively, the procedure produced the following estimates for keyblocks per square metre (with “height of exposure” = 2 and “length of exposure” = 10): 0.82, 0.78, 0.68, 0.91, 0.84, 0.71, 0.78, 0.67, 0.67, 0.63 (av. 0.75). Similarly for “mean” = 1.5 and “STD” = 0.5 the estimates obtained were 0.49, 0.36, 0.38, 0.36, 0.36, 0.34, 0.41, 0.42, 0.37, 0.41, 0.38 (av. 0.39). Now in order to evaluate the distribution of keyblock volumes we proceed as follows: Letn,=(ni,,niZ,ni3), i = 1, 2, 3, beanormal tojointsof the ith family and n, = (no,, n,, no3) be a normal to the excavation face. Consider the pyramid bounded by the planes: njfX
+nfzY
+nfjZ
=O,
i = 1,2,3, and the plane
n,,x+n,,y+n,,z=l. This pyramid has a vertex at 0 = (0,0,O)and three other vertices v, = (xi, yi, z,), i = 1,2, 3, which can be calculated by solving three systems of linear equations. We can also calculate the angles yi of the pyramid adjacent to the vertex 0: 13.00m
Length of borehola Length
ot
Height
of exposure
6.23m
Height
of tunnel
2.00m
Dip of excavation Dip direction Fault
Dip
(12)
6.00m
exposure
Searing
90.00”
face
of excavation
No in borehole --36 6 :
face
No on exooswe 55 6 65
Fig. 5
110.00”
Tracelength mean STD 3.0m 1.0m 3.0m l.Om 3.0m l.Om
Probability 0.3154 0.4609 0.2010 0.2541
Base anale 129’ 23” :to 0
i=1,2,3. All considered keyblocks are similar to each other and to the pyramid constructed above. Therefore y, represents the angles of these pyramids adjacent to the vertex inside the rock. Moreover, let S be the area of the pyramid base visible on the excavation face. Then the volume V of this pyramid is: V = KS3’2,
(13)
SHAPIRO and DELPORTz Cumulative distribution
I
0
0.04
381
of keyblock volumrs
I
I
0.12
0.16
I
0.06
STATISTICAL ANALYSIS OF JOINTED ROCK DATA
I
0.20
I
0.24
Volume Fig. 6
where the coefficient K is the same for all considered keyblock pyramids and can be calculated as follows: Consider the constructed pyramid with vertices v, , v2, v3 and 0. The volume V* of this pyramid is: XI
det i
YI
ZI
x2
Y2
z2
x3
Y3
z3
1
and the height h* is 1
h* =.
_ JG,
_ + nt
_* + nb
(The height h * is equal to one if ]ln, 1)= 1.) Then: I/* = (fh*)“z(V*)‘/’
(14)
Formula (13) allows quick calculation of keyblock volumes. Volumes of the identified keyblocks, simulated by the procedure, were calculated and stored in the computer memory. The cumulative distribution of these volumes was then calculated. An example of such distribution (for the tunnel wall of the case described in Fig. 5) is shown in Fig. 6. The cumulative distribution of Fig. 6 indicates, for instance, that about 90% of the keyblocks have a volume less than 0.06 m3. In the remainder of this section the problem of estimating the probability that a suspected keyblock may form a complete keyblock is discussed. Consider a line segment AB which forms side QR of the triangle AQPR representing the base of a suspected keyblock. The line segment AB represents the trace of intersection of a joint with the excavation face. This joint forms
Fig. 8
a face QRC of the keyblock pyramid with the vertex C inside the rock. The shape of joint is idealized to form a circle and it is assumed that the plane of the excavation face intersects them at “random”. We have to calculate the probability pi that the triangle AQRC constructed on the interval QR of the “random” chord AB falls entirely within the circle. Of course, this probability depends on the distribution of the “random” chord AB and distances AQ = bi, QR = c, and RB = di. We assume that the intersection of the excavation plane and the joint is arbitrary, and that therefore the intersection D of AB with the perpendicular diameter is uniformly distributed on that diameter. Then pi is given by the proportion of the diameter cut off by the chord AB when the point C lies on the circumference of the circle. After certain calculations we obtain: pi = [sin(8,/2)12, where Ji is the angle L ACB which can be calculated as follows. Consider the angles:
A 0
Fig. I
Fig. 9
SHAIIIP.0 and DELPORT:
?82
STATJSTJCAL ANALYSIS OF JOINTED ROCK DATA
the same family. However, by definition, such a quadrangle cannot be a base of a keyblock (see [3]). Extension of the estimation procedure to keyblocks formed by families of four or more joints, is a matter for further development. Notice that once the model simulation, described in the third section, is completed it basically becomes the question of computer time to count ail simulated keyblocks formed by families of, say, four joints and to calculate the corresponding estimates.
and /I& =
L
CRQ =
7T -
a, --
3’, .
Then
where sin ii = dJ(sin a,/sin yi)‘cf + df - 2cidi(sin a,/sin y,jcos(a, + yi)]-“2 Sin(a, + ri)
Acknowledgements-This paper is based on a project supported by the Centre for Advanced Computing and Decision Support, CSIR, Pretoria and follows CSIR Technical Report PKOMP 89/2. The authors wish to express their thanks to H. W. Ittman for careful management of the project, to F. A. Vreede for helpful suggestions and to the referee for constructive comments.
and sin vi = b,[(sin /$/sin yi)2cf + bf - 2cibi(sin /Ii/Sin r,)COS(/3i+ yi)]-“’ Sin(fii + ri). If we assume that joints of different families are independently distributed, then the probability P that a suspected keyblock forms a complete keyblock is p,p2pj. It follows that under the above model: P = [sin(s,/2)sin(82/2)sin(6,/2)]2.
Notice that the angles ai and hence the probability depend on the distances bi, Ci, 4, i = 1,2,3. CONCLUDING
(15) P
REMARKS
It should be understood that the procedure developed here gives only an approximation of reality. Large errors can occur due to inaccurate estimates of the parameters p and r~ of the log-normal distribution of trace lengths, because for instance, line segments of an identified family are not exactly parallel. It is recommended that the program should be run several times with slightly different values of the input parameters, in order to get some idea about variability of the output estimators. It is emphasized again that all suspected keyblocks visible on the simulated excavation face should be counted even if they overlap each other (which is often the case). Notice that it is quite possible for line segments from three families under consideration to form a closed quadrangle with two sides formed by parallel lines from
Accepted for publication 25 January 1991.
REFERENCES 1. Shi G. H. and Goodman R. E. Underground support design using block theory to determine keyblock bolting requirements. Proc. SANGORM
Symp. on Rock Mech.
in the Design of Tunnels,
Pretoria, pp. 81-105 (1963). 2. Shi G. H. and Goodman R. E. Keyblock bolting. Proc. Inr. Symp. on Rock Bolting, Abisko, pp. 143-167 (1983). 3. Delport J. L. and Martin D. H. A multiplier method for identifying keyblocks in excavations through jointed rock. SIAM J. Aig. Disc. Meth. 7, 321-330 (1986). 4. Vreede F. A. and Martin D. H. Keyblock characterization of jointed rock. SANGORM Symp. on Rock Mass Characterization (1985).
5. McCullagh P. and Lang P. Stochastic models for rock instability in tunnels. J. R. Statist. Sot. B. 46, 344-352 (1984). 6. Priest S. D. and Hudson J. A. Estimation of discontinuity spacing and trace length using scanline surveys. ht. J. Rock. Mech. Mm. Sci. & Geomech. Abstr. 18, 183-197 (1981).
7. Pahl P. J. Estimating the mean length of discontinuity traces. Inf. J. Rock Mech. Min. Sci. & Geomech. Abstr. X3, 221-228 (1981). 8. Rouleau A. and Gale J. E. Statistical characterization of the fracture system in the Stripa granite. Sweden. Int. J. Rock Mech. Min. Sci. & Geomech. Abstr. 22, 353-367 (1985). 9. Bury K. V. Statistical Models in Applied Sciences. Wiley, New York (1975). 10. Kendall M. G. and Moran P. A. P. Geometrical Probability. Griffin’s Statistical Monographs, London (1963). 11. Galambos J. The Asymptotic Theory of Extreme Order Statistics. Wiley, New York (1978).