gmGeostats (version 0.11.3)

anis_GSLIBpar2A: Produce anisotropy scaling matrix from angle and anisotropy ratios

Description

Produce anisotropy matrix (as the transposed of the Cholesky decomposition) from angle and anisotropy ratios

Usage

anis_GSLIBpar2A(ratios, angles, inv = FALSE)
anis2D_par2A(ratio, angle, inv = FALSE)
anis3D_par2A(ratios, angles, inv = FALSE)

Value

a 3x3 matrix of anisotropy.

If inv=TRUE (or 1) the output is a matrix A such that norm(h %*% A)

allows to use isotropic variograms, being h = c(hx, hy, hz) the lag vectors.

If inv=FALSE (or 0) the output is a matrix A such that norm(h %*% solve(A))

allows to use isotropic variograms.

Other values are meaningless.

Arguments

ratios

vector of two values between 0 and 1 giving the anisotropy ratios of medium/largest smallest/largest ranges

angles

as defined in gstat::vgm (and indeed GSLIB). For anis2D_par2A 'angle' is the direction of maximum range, i.e. largest spatial continuity, measured clockwise from North

inv

boolean or integer, see return for details

ratio

an anisotropy ratio (min/max range)

angle

direction of maximum range, i.e. largest spatial continuity, measured clockwise from North

Functions

• anis2D_par2A: 2D case

• anis3D_par2A: 3D case

See Also

Other anisotropy: AnisotropyRangeMatrix(), AnisotropyScaling(), as.AnisotropyRangeMatrix(), as.AnisotropyScaling(), is.anisotropySpecification()

Examples

## ratio=0.5, azimuth 30?? (i.e. direction 60??)
A = anis2D_par2A(1, 30)
A
AAt = A %*% t(A)
 #  project the bisector 1:1 (i.e. 45??)
(k = c(1,1,0) %*% A)
atan2(k[2], k[1]) * 180/pi  # should be 15
sqrt(sum(k^2))
sqrt( c(1,1,0) %*% AAt %*% c(1,1,0) )
A = anis2D_par2A(0.5, 60)
rd = 60 * pi/180
A
A %*% t(A)
c(cos(rd), sin(rd),0) %*% A #  should be 1
c(-sin(rd), cos(rd),0) %*% A #  should be +/- sqrt(2)
c60 = cos(60*pi/180)
s60 = sin(60*pi/180)
c30 = cos(30*pi/180)
s30 = sin(30*pi/180)
#  in the new coordinates, 60cwN is (0,1,0)
R60p = anis3D_par2A(ratios=c(1,1), angles=c(60,0,0))
c(s60, c60, 0) %*% R60p
R6030 = anis3D_par2A(ratios=c(1,1), angles=c(60,30,0))
# the original X axis is positive on newX and newY, but negative on newZ
c(1,0,0) %*% R6030
# rotate first direction 60 degrees azimuth, then dip 30degrees upwards
c( c60*c30, -s60*c30, s30) %*% R6030
(Ranis = anis3D_par2A(ratios=c(0.5,0.25), angles=c(60,30,0)) )