subroutine aptrklc (au, av, abu, abv, cu, cv, du, dv, & dintmn, dintmx, np, tol, & nint, eu, ev, dint, nerr) ccbeg. cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc c c SUBROUTINE APTRKLC c c call aptrklc (au, av, abu, abv, cu, cv, du, dv, c & dintmn, dintmx, np, tol, c & nint, eu, ev, dint, nerr) c c Version: aptrklc Updated 1990 November 29 17:40. c aptrklc Originated 1990 January 11 15:10. c c Author: Arthur L. Edwards, LLNL, L-298, Telephone (925) 422-4123. c c c Purpose: To find, for each of the np sets of input data, any c intersection point e = (eu, ev) of the linear track through c point a = (au, av) with the direction vector ab = (abu, abv), c and the line segment from point c = (cu, cv) to point c d = (du, dv), all in the uv plane of uvw space, for which c (1) the distance from point "a" to point "e" is between the c limits dintmn and dintmx, (2) point "e" is between point "c" c and "d", and (3) the track crosses line "cd" from left to right c in the uv plane. This is a zone exit if the points "c" and "d" c are the vertices of a zone edge, counterclockwise around the c zone in the uv plane. Flag nint indicates the type of c intersection found. The distance dint of the intersection c from point "a" is also returned. c Flag nerr indicates any input error. c c History: 1990 March 14. Changed tol to 0.0 in call to unit vector c subroutine. Allows small magnitudes. c c Input: au, av, abu, abv, cu, cv, du, dv, dintmn, dintmx, np, tol. c c Output: nint, eu, ev, dint, nerr. c c Calls: aptfdav, aptvadc, aptvaxc, aptvdic, aptvubc c c c Glossary: c c abu, abv Input The u and v components of 2-D direction vector "ab". c Must not both be zero. Size np. c c au, av Input The u and v coordinates of point "a". Size np. c In the uv plane. c c cu, cv Input The u and v coordinates of point "c" in the uv plane. c Must differ from "d", based on tol. Size np. c c dint Output The distance of the point of intersection "e" from c point "a". Positive if in the same direction as c vector "ab". Size np. c c dintmn Input Minimum allowed value of distance to intersection. c Size np. c c dintmx Input Maximum allowed value of distance to intersection. c Size np. c c du, dv Input The u and v coordinates of point "d" in the uv plane. c Must differ from "c", based on tol. Size np. c c eu, ev Output The u and v coordinates of the point of intersection c of the line through point "a" with direction c vector "ab", and line "cd", all in the uv plane. c Size np. c c nerr Output Indicates an input error, if not 0. c 1 if np is not positive. c c nint Output Number of acceptable intersections found. Size np. c 0 if none. c 1 if an acceptable intersection was found. c 2 if the track coincides with part of line "cd". c c np Input Size of arrays. 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---- Component u of vector "ab". dimension abu (1) c---- Component v of vector "ab". dimension abv (1) c---- Coordinate u of point "a". dimension au (1) c---- Coordinate v of point "a". dimension av (1) c---- Coordinate u of point "c". dimension cu (1) c---- Coordinate v of point "c". dimension cv (1) c---- Distance from point "a" to point "e". dimension dint (1) c---- Minimum distance to intersection. dimension dintmn (1) c---- Maximum distance to intersection. dimension dintmx (1) c---- Coordinate u of point "d". dimension du (1) c---- Coordinate v of point "d". dimension dv (1) c---- Coordinate u of point "e". dimension eu (1) c---- Coordinate v of point "e". dimension ev (1) c---- Number of acceptable intersections. dimension nint (1) c.... Local variables. c---- Component u of vector "ac". common /laptrklc/ acu (64) c---- Component v of vector "ac". common /laptrklc/ acv (64) c---- Distance from point "a" to line "cd". common /laptrklc/ atocd (64) c---- Distance from point "c" to line "ab". common /laptrklc/ ctoab (64) c---- Component u of vector "cd". common /laptrklc/ cdu (64) c---- Component v of vector "cd". common /laptrklc/ cdv (64) c---- Fractional distance of "e" along "cd". common /laptrklc/ fint (64) c---- A very small number. common /laptrklc/ fuz c---- Component u of unit vector "hab". common /laptrklc/ habu (64) c---- Component v of unit vector "hab". common /laptrklc/ habv (64) c---- Component u of unit vector "hcd". common /laptrklc/ hcdu (64) c---- Component v of unit vector "hcd". common /laptrklc/ hcdv (64) c---- Index in arrays. common /laptrklc/ n c---- First index of subset of data. common /laptrklc/ n1 c---- Last index of subset of data. common /laptrklc/ n2 c---- Limit flag from aptfadv. common /laptrklc/ nlim c---- Index in external array. common /laptrklc/ nn c---- Size of current subset of data. common /laptrklc/ ns c---- Sine of angle between "ab" and "cd". common /laptrklc/ sinth (64) c---- Magnitude of a vector. common /laptrklc/ vlen (64) c---- Magnitude of vector "ab". common /laptrklc/ vlenab (64) c---- Magnitude of vector "ac". common /laptrklc/ vlenac (64) c---- Magnitude of vector "cd". common /laptrklc/ vlencd (64) cbugc***DEBUG begins. cbug 9901 format (/ 'aptrklc finding intersection of the line through:' / cbug & (i3,' au,av= ',1p2e22.14 / cbug & ' with direction vector:' / cbug & ' abu,abv=',1p2e22.14 / cbug & ' with the line through:' / cbug & ' cu,cv= ',1p2e22.14 / cbug & ' du,dv= ',1p2e22.14 / cbug & ' dintmn= ',1pe22.14,' dintmx=',1pe22.14)) cbug write ( 3, 9901) (n, au(n), av(n), abu(n), abv(n), cbug & cu(n), cv(n), du(n), dv(n), dintmn(n), dintmx(n), cbug & n = 1, np) cbugc***DEBUG ends. c.... Initialize. c---- A very small number. fuz = 1.e-99 nerr = 0 c.... Test for input errors. if (np .le. 0) then nerr = 1 go to 210 endif c.... Set up the indices of the first subset of data. n1 = 1 n2 = min (np, 64) 110 ns = n2 - n1 + 1 c.... Find the unit vector "hab" parallel to 2-D direction vector "ab". call aptvubc (abu(n1), abv(n1), ns, 0., & habu, habv, vlenab, nerr) c.... Find the 2-D vector "cd" from point "c" to point "d". call aptvdic (cu(n1), cv(n1), du(n1), dv(n1), ns, tol, & cdu, cdv, vlencd, nerr) c.... Find the unit vector "hcd" parallel to 2-D direction vector "cd". call aptvubc (cdu, cdv, ns, 0., & hcdu, hcdv, vlencd, nerr) c.... Find the 2-D vector "ac" from point "a" to point "c". call aptvdic (au(n1), av(n1), cu(n1), cv(n1), ns, & tol, acu, acv, vlenac, nerr) c.... Find the perpendicular distance from point "a" to line "cd". call aptvaxc (acu, acv, hcdu, hcdv, ns, tol, & atocd, nerr) c.... Find the perpendicular distance from point "c" to line "ab". call aptvaxc (acu, acv, habu, habv, ns, tol, & ctoab, nerr) c.... Find the component of unit vector "hab" perpendicular to line "cd". call aptvaxc (habu, habv, hcdu, hcdv, ns, tol, sinth, nerr) c.... Find the distance of the intersection point "e" from point "a", c.... and the fractional distance of point "e" along line "cd". c---- Loop over subset of data. do 120 n = 1, ns nn = n + n1 - 1 dint(nn) = atocd(n) / (sinth(n) + fuz) fint(n) = ctoab(n) / (sinth(n) * vlencd(n) + fuz) c---- End of loop over subset of data. 120 continue c.... Impose limits on the values of fint. call aptfdav (fint, ns, 1, tol, nlim, nerr) c.... Find the the point of intersection "e". call aptvadc (cu(n1), cv(n1), 1., fint, cdu, cdv, ns, tol, & eu(n1), ev(n1), vlen, nerr) c.... Test the intersection for acceptability. c---- Loop over subset of data. do 130 n = 1, ns nn = n + n1 - 1 if ((dint(nn) .ge. dintmn(nn)) .and. & (dint(nn) .le. dintmx(nn)) .and. & (fint(n) .ge. 0.0) .and. & (fint(n) .le. 1.0) .and. & (sinth(n) .gt. 0.0) ) then nint(nn) = 1 else nint(nn) = 0 endif c---- End of loop over subset of data. 130 continue c.... See if all data subsets are done. c---- Do another subset of data. if (n2 .lt. np) then n1 = n2 + 1 n2 = min (np, n1 + 63) go to 110 endif cbugc***DEBUG begins. cbug 9902 format (/ 'aptrklc results:' / cbug & (i3,' nint=',i2 / cbug & ' dint= ',1pe22.14 / cbug & ' eu,ev= ',1p2e22.14)) cbug write ( 3, 9902) (n, nint(n), dint(n), eu(n), ev(n), n = 1, np) cbugc***DEBUG ends. 210 return c.... End of subroutine aptrklc. (+1 line.) end UCRL-WEB-209832