subroutine aptplsp (rsph, px, py, pz, ax, ay, az, bx, by, bz, np,
& tol, dpmin, rcir, cx, cy, cz, itrun, nerr)
ccbeg.
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
c
c SUBROUTINE APTPLSP
c
c call aptplsp (rsph, px, py, pz, ax, ay, az, bx, by, bz, np,
c & tol, dpmin, rcir, cx, cy, cz, itrun, nerr)
c
c Version: aptplsp Updated 1990 November 28 14:50.
c aptplsp Originated 1990 March 16 14:10.
c
c Author: Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123.
c
c
c Purpose: To find, for each of np sets of input data, the minimum distance
c dpmin to the point p = (px, py, pz), from the plane through
c the point a = (ax, ay, az) with normal vector b = (bx, by, bz),
c and the coordinates c = (cx, cy, cz) of the point in the plane
c nearest point "p". In addition, to find the radius rcir of any
c circle of intersection of the plane with the sphere centered
c at point "p" with radius rsph.
c Flag itrun indicates truncation of dpmin to zero (1) or
c too small a magnitude of vector "b" (2).
c Flag nerr indicates any input error.
c
c Input: rsph, px, py, pz, ax, ay, az, bx, by, bz, np, tol.
c
c Output: dpmin, rcir, cx, cy, cz, itrun, nerr.
c
c Calls: aptptpl
c
c Glossary:
c
c ax,ay,az Input The x, y, z coordinates of point "a" in the plane.
c Size np.
c
c bx,by,bz Input The x, y, z components of vector "b" normal to the
c plane. Magnitude must exceed tol. Size np.
c
c cx,cy,cz Output The x, y, z coordinates of the point in the plane
c nearest point "p". Size np.
c The center of any circle of intersection.
c Returned as point "p" if normal vector "b" is too
c short, based on tol (itrun = 2).
c
c dpmin Output The perpendicular distance to point "p" from the
c plane through point "a" with normal vector "b".
c Positive if point "p" is in the same direction
c from the plane as the normal vector "b".
c Truncated to zero if less than the estimated error
c in its calculation, based on tol (itrun = 1).
c Returned as zero if normal vector "b" is too short,
c based on tol (itrun = 2). Size np.
c
c itrun Output Indicates a special result for one data set:
c 0 if dpmin is not zero, and no input errors occur.
c 1 if the value of dpmin is truncated to zero, when
c less than the estimated error in its calculation,
c based on tol.
c 2 if normal vector "b" is too short, based on tol.
c The orientation of the plane is unknown, and
c dpmin and point "c" cannot be calculated.
c Size np.
c
c nerr Output Indicates an input error, if not 0.
c 1 if np is not positive.
c
c np Input Size of arrays.
c
c px,py,pz Input The x, y, z coordinates of point "p". Size np.
c The center of the sphere with radius rsph.
c
c rcir Output The radius of any circle of intersection centered at
c point "c". Size np. Positive if intersection
c occurs, zero if tangency occurs, negative if no
c intersection or tangency occurs, or if itrun = 2
c or 3.
c
c rsph Input The radius of the sphere centered at point "p".
c Size np. The absolute value is used.
c
c tol Input Numerical tolerance limit.
c On Cray computers, recommend 1.e-5 to 1.e-11.
c
cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
ccend.
c.... Dimensioned arguments.
c---- Coordinate x of point "a".
dimension ax (1)
c---- Coordinate y of point "a".
dimension ay (1)
c---- Coordinate z of point "a".
dimension az (1)
c---- Component x of vector "b".
dimension bx (1)
c---- Component y of vector "b".
dimension by (1)
c---- Component z of vector "b".
dimension bz (1)
c---- Coordinate x of point "c".
dimension cx (1)
c---- Coordinate y of point "c".
dimension cy (1)
c---- Coordinate z of point "c".
dimension cz (1)
c---- Distance from plane to point "p".
dimension dpmin (1)
c---- 1 if dpmin truncated to zero.
dimension itrun (1)
c---- Coordinate x of point "p".
dimension px (1)
c---- Coordinate y of point "p".
dimension py (1)
c---- Coordinate z of point "p".
dimension pz (1)
c---- Radius of circle at "c".
dimension rcir (1)
c---- Radius of sphere at "p".
dimension rsph (1)
c.... Local variables.
c---- Equals abs(rsph) - abs(dpmin).
common /laptplsp/ dr
c---- Equals abs(rsph) + abs(dpmin).
common /laptplsp/ ds
c---- Index in arrays.
common /laptplsp/ n
cbugc***DEBUG begins.
cbug 9901 format (/ 'aptplsp finding intersection of the sphere at:')
cbug 9902 format (i3,' px,py,pz=',1p3e22.14 /
cbug & ' rsph= ',1pe22.14 /
cbug & ' with the plane through:' /
cbug & ' ax,ay,az=',1p3e22.14 /
cbug & ' with normal vector:' /
cbug & ' bx,by,bz=',1p3e22.14)
cbug write ( 3, 9901)
cbug write ( 3, 9902) (n, px(n), py(n), pz(n), rsph(n),
cbug & ax(n), ay(n), az(n),
cbug & bx(n), by(n), bz(n), n = 1, np)
cbugc***DEBUG ends.
c.... Initialize.
nerr = 0
c.... Test for input errors.
if (np .le. 0) then
nerr = 1
go to 210
endif
c.... Find the distance dpmin and coordinates "c" of the point in
c.... the plane nearest point "p".
call aptptpl (px, py, pz, ax, ay, az, bx, by, bz, np, tol,
& dpmin, cx, cy, cz, itrun, nerr)
c.... Find the radius of the circle of intersection.
c---- Loop over subset of data.
do 120 n = 1, np
dr = abs (rsph(n)) - abs (dpmin(n))
ds = abs (rsph(n)) + abs (dpmin(n))
if (abs (dr) .lt. tol * ds) then
dr = 0.0
endif
rcir(n) = sign (sqrt (abs (dr * ds)), dr)
if (itrun(n) .gt. 1) then
rcir(n) = -abs (rcir(n))
endif
c---- End of loop over subset of data.
120 continue
cbugc***DEBUG begins.
cbug 9903 format (/ 'aptplsp results:' /
cbug & (i3,' dpmin= ',1pe22.14,' itrun=',i2 /
cbug & ' rcir= ',1pe22.14 /
cbug & ' cx,cy,cz=',1p3e22.14))
cbug write ( 3, 9903) (n, dpmin(n), itrun(n), rcir(n),
cbug & cx(n), cy(n), cz(n), n = 1, np)
cbugc***DEBUG ends.
210 return
c.... End of subroutine aptplsp. (+1 line.)
end
UCRL-WEB-209832